Françoise Olivier Hyrien

Benoit Le Tallec

Magali Hennion

Elizaveta Laurent Lacroix

Jean-Michel Barbara

Jérémy Marion

Vincent Yacouba

Marcelo Laurence

Eric, Salambo Roberto Carlos

Soha, Sarah, Haydar, Zahraa, Riva, Batoul, Zeinab Abed Orsola

Hussein El Salam, Hiba

Keywords: ascending classication, Correlation matrices representation and comparison

aurait jamais pu être réalisé.

Résumé

Le programme spatiotemporel de réplication de l'ADN est régulé au cours du développement et altéré durant la progression cancéreuse. Nous proposons une caractérisation originale de la plasticité du programme de réplication de l'ADN se basant sur le prolage de 12 lignées cellulaires humaines normales ou cancéreuses par la méthode OK-seq de purication et séquençage des fragments d'Okazaki qui permet de déterminer l'orientation de la progression des fourches de réplication (OFR) à haute résolution (10 kilo bases). L'analyse comparative des prols OFR montre que les changements réplicatifs permettent la classication des lignées cellulaires en fonction de leur tissu d'origine, la nature cancéreuse ou non de la lignée n'intervenant qu'en second ordre. Il n'apparait pas de point chaud pour l'accumulation des changements réplicatifs, ceux-ci étant largement dispersés sur tout le génome. Néanmoins, les régions riches en G+C et en gènes actifs, répliquées précocement au cours de la phase S, ont le programme de réplication le plus stable, elles présentent une forte densité de zones d'initiation de la réplication (ZI) ecaces et conservées entre lignées cellulaires. En contraste, les dernières régions répliquées, à faible densité de gènes et pauvres en G+C, présentent peu de ZI ecaces, souvent spéciques d'un tissu ou d'une lignée. Ceci nous conduit à quantier le degré de dissociation entre ZI et activation de la transcription. Ce travail propose un panorama original des modications du programme de réplication au cours de la diérentiation normale ou pathologique, dont un contrôle lignée cellulaire spécique des ZI dans les déserts de gènes à réplication tardive.

Abstract

The spatiotemporal program of DNA replication is regulated during development and altered during cancer progression. We propose an original characterization of the plasticity of the DNA replication program based on the proling of 12 normal or cancerous human cell lines by the OK-seq method of purication and sequencing of Okazaki fragments which allows to determine the orientation of the progression of replication forks (RFD) at high resolution (10 kilo bases). Comparative analysis of the RFD proles shows that the replicative changes allow the classication of the cell lines according to their tissue of origin, the cancerous or noncancerous nature of the cell line intervening only in second order. There is no hotspot for the accumulation of replicative changes, they are widely dispersed throughout the genome.

Nevertheless, the G+C rich and active gene regions, replicated early in the S phase, have the most stable replication program, they present a high density of ecient replication initiation zones (IZ) conserved between cell lines. In contrast, the late replicated, low gene density and low G+C content regions have few ecient IZs, often specic to a tissue or lineage. This leads us to quantify the degree of dissociation between IZ and activation of transcription. This work provides an original overview of replication program changes during normal or pathological dierentiation, including a cell line specic control of IZ in late-replication gene deserts.

Introduction

Life depends on the cells ability to preserve, retrieve and translate the genetic instructions needed to build a living organism and keep it alive. This hereditary information is transmitted from a cell to its daughter cells at the time of cell division and, for one organism, from one generation to another by the reproductive cells of this organism. In each living cell, these instructions are carried by the genes, the elements containing the information which determines the characteristics of a species as well as those specic to each individual.

The biological carriers of the genetic information are deoxyribonucleic acid (DNA) molecules.

Hence, DNA replication, the duplication of these molecules, is an essential mechanism for the living beings. This mechanism must be as faithful as possible, so that the genetic heritage that will be transmitted to the daughter cells is almost perfectly identical to that of the mother cell.

Indeed, the short-term survival of an organism depends on the prevention of modication of its genetic content. DNA modications such as mutations can be extremely deleterious events, associated with pathologies like cancer. However, there may be benec mutations responsible for the long term evolution of the species. During the S phases of cell cycle where DNA replication occurs in eukaryotes, the genome is particularly vulnerable and many types of lesions can block the progression of the replication machinery. In fact, prolonged stopping of the replication machinery can lead to its collapse, double-strand breaks, and genetic instability.

Stabilization and backup of blocked machinery and nally the completion of replication are essential for cell survival and genome preservation [1,2].

In this context, the present work aim was to contribute to the better understanding of DNA replication process in human by identifying and characterizing its modications observed between 12 human cell lines. This work is based on the recent development of a proling protocol allowing to measure genome-wide the so-called Replication Fork Directionality (RFD) [3]. In this introductory chapter, we will discuss how the discovery of the structure of DNA led to the modeling of the mechanism of DNA replication, and the current experimental methods to study DNA replication. Finally, we will describe the association between the organization of DNA replication and other genomic features.

1.1 Deoxyribonucleic acid structure and copying mechanism

In the 1940s, there has been a vivid devate among biologists to recognize DNA as the preincipal carrier of genetic information. This long polymer macromolecule was considered as too simple, consisting only of four types of chemically similar nucleotides. These nucleotides are composed of a phosphate group (or phosphoric acid), an aldopentose (2'-deoxy-D-ribose) and a nitrogen heterocyclic base, either one of two pyrimidines: cytosine (C) and thymine (T), or two purines: adenine (A) and guanine (G) [4]. In the early 1950s, for the rst time, DNA was analyzed by X-ray diraction, a technique that makes it possible to determine the threedimensional atomic structure of a molecule. The rst results of the X-ray diraction showed that the DNA was very likely constituted of two co-axial molecules wound into a helix [5]. The fact that DNA is double-stranded constitutes one of the main hypothesis that led to the model of DNA structure published in 1953 by Watson and Crick [START_REF] Watson | Molecular structure of nucleic acids[END_REF]. In this paper, they described a structure based on specic base pairing of adenine with thymine and guanine with cytosine at the center of a double helix and holding together two single strands of DNA. The sugarphosphate backbone of each strand of a DNA double helix has a unique chemical direction, or polarity, determined by the manner in which each sugar residue is bonded to the next, and the two strands of the double helix are anti-parallel [START_REF] Fangman | Activation of replication origins within yeast chromosomes[END_REF]; that is, they go in opposite directions.

Specic base complementarity between the two DNA stands implies that each strand contains all the information to reconstruct the whole double-stranded molecule. At the very end of this brief scientic blockbuster, Watson and Crick comment almost in an aside: It has not escaped our notice that the specic pairing we have postulated immediately suggests a possible copying mechanism for the genetic material . The structure of DNA would allow for its exact duplication, an essential property required for the genetic material.

One month after their historical article in Nature, Watson and Crick published a second article [START_REF] Watson | Genetical implications of the structure of deoxyribonucleic acid[END_REF], suggesting how DNA could be duplicated. In this article, they proposed that each strand serves as a model for the synthesis of a complementary daughter strand. For this duplication to occur, the two parental strands must separate, requiring a rupture of the hydrogen bonds between paired bases, and unfolding of the paired strands. In this model, called semiconservative replication, each new molecule of DNA is composed of a strand derived from the original parent molecule and a newly synthesized strand (Figure 1.1). The discovery of the structure of DNA immediately suggested responses on how the information needed to dene an organism can be duplicated and copied from generation to generation.

The semiconservative model indeed describes the main steps of DNA replication in vivo. DNA replication, also known as DNA duplication, is the process in which DNA is synthesized in vivo. This mechanism makes it possible to obtain, from one DNA molecule, two molecules identical to the initial molecule, for distribution to the two daughter cells during cell division. Because nucleotide A will successfully pair only with T, and G with C, each strand of a DNA double helix labeled here as the S strand and its complementary S' strand can serve as a template to specify the sequence of nucleotides in its complementary strand. In this way, both strands of a DNA double helix can be copied precisely. In this model, the two new molecules are formed by one original (template) strand and a newly synthesized strand, this model is termed semiconservative. (Adapted from Essential Cell Biology fourth edition [START_REF] Alberts | Essential cell biology[END_REF]).

However, it only captures part of the complexity involved in regulating the faithful replication of a complete genome, like the human genome made of 6.6 billion nucleotides distributed on 23 pairs of chromosomes, which requires further levels of complexity. According to the so-called "replicon" model proposed in 1963 by François Jacob, François Cuzin and Sydney Brenner [START_REF] Jacob | On the regulation of DNA replication in bacteria[END_REF], replication of a chromosome starts with the binding of some "initiator" protein to a specic "replicator" DNA sequence called origin of replication. The replicon model is the paradigm of the replication of circular, prokaryotic 1 chromosomes having a single replication origin (Figure 1.2). The initiator-replicator association makes it possible to open the DNA locally and trigger replication. The recruitment of additional factors initiates the bi-directional progression along the chromosome of two divergent replication forks each assuring the simultaneous polymerization of the two new DNA strands as predicted by the semiconservative replication model (Figure 1.2). This forms a replication bubble until the two meets at the termination site, opposite from the origin on the chromosome (Figure 1.2). 1 A prokaryote is a living being whose cell structure does not contain a nucleus, and almost never any membranous organelles Figure 1.2: Semiconservative model of bacterial DNA replication. The circular double-stranded DNA chromosome of most bacteria is duplicated by bidirectional synthesis. It begins near a special site on the chromosome, called the origin (O) where the DNA opens up as hydrogen bonds between the bases in the complementary strands are enzymically broken and the strands separate, forming a replication 'bubble' between two DNA forks as the two strands separate on either side. The DNA continues unzipping at these forks, in both directions simultaneously. These so-called replication forks move around the chromosome in opposite directions as the DNA is unzipped and duplicated (hence bidirectional synthesis). DNA polymerase reads the base sequence of the unzipped single strands of parental template DNA (blue lines) and matches the sequence of bases with new bases which are polymeried to create the complementary daughter strands (red dashed lines). This continues until the replication forks reach the terminus (T) at the far-side of the chromosome, opposite the origin. (Adapted from http://cronodon.com/ BioTech/Bacteria_Growth.html).

The situation in eukaryotes 2 bears some resemblance with prokaryotes, we have the same process of replication, but there are several origins of replication on each chromosomes as we see in Figure 1.3. Replication of a chromosome is initiated at a number of replication origins, diverging replication forks propagate until two converging forks collide at a terminus of replication resulting in replication bubble mergers until chromosome complete replication (Figure 1.3) [START_REF] Bell | DNA replication in eukaryotic cells[END_REF][START_REF] Depamphilis | Regulating the licensing of DNA replication origins in metazoa[END_REF].

One of the specicity of the Eukaryotic cells is the cell cycle during which two phases are repeated cyclically [START_REF] Depamphilis | Genome duplication (concepts, mechanisms, evolution, and disease) garland science[END_REF][START_REF] Diey | The chromosome replication cycle[END_REF]. Mitosis corresponding to the division of a mother cell into two daughter cells strictly identical genetically, and the interphase which is subdivided in 3 subphases: G1, S and G2. In G1 phase, the cell performs its normal metabolism and grows to a critical size triggering S phase. In the S phase, DNA replication occurs. Finally, in the G2 phase, the cellular growth is over and the centrosomes replicate, to allow the smooth course of the mitosis. For all organisms, the DNA replication during S phase is subdivided into three stages: the initiation which, through the opening of the DNA at the level of the origins, allows the recruitment of the polymerization machinery and the start of the DNA synthesis, the elongation which corresponds to the continuation of the synthesis, and the termination which consists in the fusion of two convergent forks. 3) after the time of the micrograph and after the merger of the two rightmost replication bubbles despicted in (2). The orange lines represent the two parental DNA strands; the red lines represent the newly synthesized DNA strands. (Adapted from Essential Cell Biology fourth edition [START_REF] Alberts | Essential cell biology[END_REF]).

.

Replication initiation

Experiments carried out using unicellular eukaryotic models, such as the yeasts Saccharomyces cerevisiae and Schizosaccharomyces pombe, made it possible to identify the protein actors that act at the level of the origins of replication to initiate the synthesis of DNA [1416]. The activity of these initiation proteins is nely regulated by specic proteins, such as kinases, to ensure a unique duplication of the genome at each cell cycle and thus guarantees its stability.

The goal of the initiation of replication is to set up a large number of factors, proteins necessary for the elongation of replication. It includes during late M and early G1 phases loading on the DNA of the ORC 3 (Origin Recognition Complex) protein complex and the assembly of the pre-replication replication forks and between 30 000 to 50 000 origins [START_REF] Méchali | Eukaryotic DNA replication origins: many choices for appropriate answers[END_REF][START_REF] Cvetic | Eukaryotic origins of DNA replication: could you please be more specic?[END_REF], because timely duplication of large linear chromosomes requires establishment of replication forks at multiple locations to complete replication in a typical 8 hours S-phase in human.

Elongation asymmetry

To complete the replication there is an essential step which is the propagation of the replication forks named elongation. Upon initiation of a replication origin, the two forks move away from the origin in opposite directions, unfolding the DNA double helix and replicating the DNA as it progresses (Figure 1.4). Forks move very fast in bacteria, about 1000 pairs of nucleotides per second, and about 10 times slower in eukaryotes. The slowness of fork movement in eukaryotes may be due to the diculties of DNA replication through the more complex chromatin structure of eukaryotic chromosomes. All DNA polymerases add new nucleotides only to the 3' end of a DNA strand. As a result, a new DNA chain can only be synthesized in the 5' to 3' direction. The two template strands being antiparallel, this can easily be accommodated for the synthesis based on the DNA template strand in the 3' to 5' relative to the replication fork movement. However, what happens on the opposite DNA template strand which seems to require synthesis in the reverse orientation than the fork movement? The new DNA strand which seems to grow in the wrong direction from 3' to 5' is in fact synthesized discontinuously, in small successive separated pieces. DNA polymerase moves backward with respect to the direction of replication-fork movement so that each new DNA fragment can be polymerized in the 5' to 3' direction (Figure 1.4). The resulting small DNA fragments -called Okazaki fragments -are then joined to form a new continuous strand (Figure 1.4). The strand of DNA that is synthesized discontinuously is called the lagging strand , because the backstitching imparts a slight delay to its synthesis; the complementary strand, which is synthesized continuously, is called the leading strand (Figure 1.4). This process requires fundamental actors including replicative polymerases. Polymerase δ catalyzes the extension of the continuous strand and polymerase assures the elongation of the discontinuous strand during replication [2023]. Many other factors are involved in the progression of the replication forks to ensure the proper progression of the S phase [2428]. For example, the proling cell nuclear antigen (PCNA) provides a strong association with DNA contributing to the processivity of the replicative complex (replisome) [START_REF] Mondol | PCNA accelerates the nucleotide incorporation rate by DNA polymerase δ[END_REF][START_REF] Maga | Proliferating cell nuclear antigen (PCNA): a dancer with many partners[END_REF].

Replication termination

The termination of replication occurs when two convergent forks from two active origins meet at the end of the elongation step. DNA ligase I catalyzes the formation of the phosphodiester bond between the two neo-synthesized strands of DNA. The topological constraints generated Introduction posed by the existence of such sequences is that they are more prone to "hot spots" of DNA breaks and chromosomal rearrangements [START_REF] Lambert | Checkpoint responses to replication fork barriers[END_REF].

Eukaryotic DNA replication program

It was early recognized thanks to autoradiography experiments that higher eukaryotes chromosome replication proceeds in units where DNA synthesis occur at dened time in S phase (Reviewed by Hand in 1978 [START_REF] Hand | Eucaryotic DNA: organization of the genome for replication[END_REF]). Immunouorescence staining of the synthesized DNA at the same time of S phase over two consecutive generations in mouse broblasts showed that replication proles are conserved from one generation to the next, and that the same sites are replicated at the same time [START_REF] Ma | Spatial and temporal dynamics of DNA replication sites in mammalian cells[END_REF]. These experiments corroborated that each DNA sequence is replicated at a specic time of the S phase, so that the activation time of an origin of replication is not random but conditioned by the replication time of its replication unit. More recently, methods were developed to map the mean replication timing (MRT) proles in a population of growing cells [START_REF] Farkash-Amar | Genome-wide analysis of the replication program in mammals[END_REF][START_REF] Gilbert | Evaluating genome-scale approaches to eukaryotic DNA replication[END_REF]. These direct measurements conrmed that eukaryotic chromosomes replicate in a temporal order known as the replication-timing program. In mammals, they have revealed the presence of 100 kb to Mb-scale regions with similar timing, called constant timing regions (CTRs), replicating either early or late in the S-phase by coordinated activation of multiple origins. These early and late MRT plateaus were shown to be separated by rather steep timing transition regions (TTRs) of size ranging from 0.1 to 0.6Mb and initialy presumed to be replicated unidirectionally by a single fork coming from the early domain [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF][START_REF] Farkash-Amar | Global organization of replication time zones of the mouse genome[END_REF]. It further revealed that, even though CTR locations are often conserved between cell types, MRT is cell-type-specic in mammals with up to 50% of the genome switching replication timing during development [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF][START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF][START_REF] White | DNA replication-timing analysis of human chromosome 22 at high resolution and dierent developmental states[END_REF][START_REF] Hansen | Sequencing newly replicated DNA reveals widespread plasticity in human replication timing[END_REF]. Similarly, comparative analysis of replication proles in 4

budding yeast species revealed species-specic replication patterns on top of an evolutionary conserved replication timing program despite chromosome rearrangements [START_REF] Müller | Conservation of replication timing reveals global and local regulation of replication origin activity[END_REF]. In human, MRT proles were shown to be linked with the genome organization (Figure 1.5a). The main genomic link being that early replicating regions correspond to regions with high gene density and high GC content whereas late replicating regions correspond to regions with low gene density and low GC content [START_REF] Rhind | DNA replication timing[END_REF][START_REF] Woodne | Replication timing of the human genome[END_REF].

More recently, it was experimentally noticed that the replication rate of TTRs is not always compatible with the unidirectional progression of one replication fork [START_REF] Guilbaud | Evidence for sequential and increasing activation of replication origins along replication timing gradients in the human genome[END_REF]. In these cases, the coherent wave of replication from the early to the late TTR borders necessarily implies a more complex replication program. This observation was conrmed genome wide by the analysis of MRT proles in seven human cell types including ES, somatic and HeLa cells [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF][START_REF] Audit | Megabase replication domains along the human genome: relation to chromatin structure and genome organisation[END_REF][START_REF] Audit | Multiscale analysis of genomewide replication timing proles using a wavelet-based signal-processing algorithm[END_REF]. This study indeed revealed that in each cell type, about half of the genome can be paved by replication U-domains where the MRT is U-shaped rather than by CTRs (Figure 1.6C). About half of the U-domains in one cell line is shared by at least another cell line (from 38.4% to 61%). The half of the genome that is covered by U-domains in fact corresponds to regions of high replication timing plasticity where replication domains may (i) reorganize Isochores with intermediate or mixed sequence features are frequently subject to replication timing regulation during dierentiation (speculatively labeled "Facultative heterochromatin") [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF][START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF]. (b) Changes in replication timing that traverse the middle of S phase accompany changes in subnuclear position and transcriptional potential. Replication early in S phase (patterns I and II) takes place within the interior euchromatic compartment, whereas replication later in S phase (patterns III, IV, and V) takes place at the nuclear periphery (pattern III), at the nucleolar periphery ( pattern III), and at internal blocks of heterochromatin (patterns IV and V) [START_REF] Berezney | Heterogeneity of eukaryotic replicons, replicon clusters, and replication foci[END_REF][START_REF] Mcnairn | Epigenomic replication: linking epigenetics to DNA replication[END_REF].(Adapted from Rhind et al. [START_REF] Rhind | DNA replication timing[END_REF]). 3)) along a 11.4 Mb long fragment of human chromosome 10 that contains 6 skew N-domains (horizontal black bars) bordered by 7 putative replication origins O1 to O7. Each dot corresponds to the skew calculated for a window of 1 kb of repeat-masked sequence. The colors correspond to intergenic (black), (+) genes (red) and (-) genes (blue). (B) Corresponding cumulative skew prole Σ obtained by cumulative addition of S-values along the sequence. (C) MRT proles from early, 0 to late, 1 for BG02 (green), K562 (red) and GM06990 (blue) cell lines. (D) Correlations between S and dMRT/dx, in BG02 (100 kb windows) along the 22 human autosomes; colors as in (A). (E) Average dMRT/dx proles (± SEM) in 663 skew N-domains after rescaling their length L to unity; colors as in (C). (Adapted from Baker et al. [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF]).

according to the so-called consolidation scenario, (ii) experience some boundary shift and (iii) emerge in a late replicating region [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF] as previously observed in the mouse genome during dierentiation [START_REF] Hiratani | Genome-wide dynamics of replication timing revealed by in vitro models of mouse embryogenesis[END_REF]. The early initiation zones at U-domain borders were described as master replication origins (MaOris) [START_REF] Audit | Megabase replication domains along the human genome: relation to chromatin structure and genome organisation[END_REF][START_REF] Hyrien | From simple bacterial and archaeal replicons to replication N/U-domains[END_REF] It was proposed that replication waves likely initiate and propagate from MaOris toward the domain center via a cascade of (non-independent) origin ring, possibly by fork-stimulated initiation, resulting in the observed parabolic (U-shaped) change of MRT across each domain [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF][START_REF] Guilbaud | Evidence for sequential and increasing activation of replication origins along replication timing gradients in the human genome[END_REF][START_REF] Audit | Multiscale analysis of genomewide replication timing proles using a wavelet-based signal-processing algorithm[END_REF][START_REF] Hyrien | From simple bacterial and archaeal replicons to replication N/U-domains[END_REF].

If we assume a nearly deterministic replication program, where at each cell cycle nearly the same set of replication origins is used, and where all replication origins re at specic times, the population average replication program reects faithfully what occurs in each cell cycle.

The time of replication of each region is then a function of its distance (assuming constant fork speed) from an active origin and the time at which the origin was activated (Figure 1.8A). In this model, some origins would be activated at the beginning of the S phase, others in the middle of the S phase or at the end of the S phase, dening early, intermediate or late replication domains, respectively. Thus, it is important to distinguish between the reproducible replication timing of broad regions of the genome and the specic replication timing of individual origin for which there is in fact little evidence in higher eukaryotes. So, independently from the biochemical details of replication initiation, it is fundamental to characterize the spatial distribution of the replication origins as well as their replication timing and eciency dened as the proportion of the cell cycle in which there are active. In this context, it was shown that budding yeast has well-dened, site-specic origins, many of which are ecient and re in as many as 90% of S phase [START_REF] Fangman | Activation of replication origins within yeast chromosomes[END_REF][START_REF] Newlon | Analysis of a circular derivative of saccharomyces cerevisiae chromosome III: a physical map and identication and location of ars elements[END_REF]. These characteristics lead to fairly homogeneous replication kinetics [START_REF] Raghuraman | Replication dynamics of the yeast genome[END_REF]. So, considering the DNA replication program in budding yeast as deterministic is rather consistent.

If in S. cerevisiae and S. pombe, it was shown that replication origins correspond to specic sites called ARS (for autonomous replication sequence) and characterized by a short sequence motif and islands rich in AT [6063], the situation is more complex in higher eukaryotes. Many methods were used to map the replication origins in human and other species [START_REF] Gilbert | Evaluating genome-scale approaches to eukaryotic DNA replication[END_REF]. Some approaches require to purify replication intermediates as in Short Nascent DNA sequencing (SNS-seq) [START_REF] Vassilev | Mapping an origin of DNA replication at a single-copy locus in exponentially proliferating mammalian cells[END_REF][START_REF] Burhans | Identication of an origin of bidirectional DNA replication in mammalian chromosomes[END_REF] and Bubble-seq [START_REF] Mesner | Bubble-chip analysis of human origin distributions demonstrates on a genomic scale signicant clustering into zones and signicant association with transcription[END_REF] other approaches map the early-ring human replication origins as in Optical Replication Mapping [START_REF] Klein | Genome-wide identication of early-ring human replication origins by optical replication mapping[END_REF] and Initiation site sequencing (Ini-seq) [START_REF] Langley | Genome-wide identication and characterisation of human DNA replication origins by initiation site sequencing (iniseq)[END_REF].

Even though the concordance between these methods has been the subject of discussions [START_REF] Langley | Genome-wide identication and characterisation of human DNA replication origins by initiation site sequencing (iniseq)[END_REF][START_REF] Hyrien | Peaks cloaked in the mist: the landscape of mammalian replication origins[END_REF], there is now increasing evidence that replication program is stochastic and that no two cell cycles use exactly the same set of replication origins and the same ring times [START_REF] Méchali | Eukaryotic DNA replication origins: many choices for appropriate answers[END_REF][START_REF] Hyrien | Peaks cloaked in the mist: the landscape of mammalian replication origins[END_REF][START_REF] Rhind | DNA replication timing: random thoughts about origin ring[END_REF]. In higher eukaryotes, the origins of DNA replication are activated in clusters as rst observed by autoradiography [START_REF] Huberman | On the mechanism of DNA replication in mammalian chromosomes[END_REF] and later by DNA combing methods [START_REF] Lebofsky | DNA replication origin interference increases the spacing between initiation events in human cells[END_REF], where simultaneous activations would happen from an excess of potential loci [7276]. Whatever the method used to identify the replication origin, the parameters that dene the choice of origins in higher eukaryotes are beginning to emerge. These parameters depend on the sequence, the topology of the DNA, the presence or the absence of nucleosome 6 and transcription [START_REF] Méchali | Eukaryotic DNA replication origins: many choices for appropriate answers[END_REF]. chromatin regions [START_REF] Karnani | Genomic study of replication initiation in human chromosomes reveals the inuence of transcription regulation and chromatin structure on origin selection[END_REF]. Studies in Drosophila and mouse showed that the CpG 7 islands are associated to replication origins and that these origins are preferentially concentrated in regions containing transcriptional promoters [START_REF] Méchali | Eukaryotic DNA replication origins: many choices for appropriate answers[END_REF][START_REF] Cayrou | Genome-scale analysis of metazoan replication origins reveals their organization in specic but exible sites dened by conserved features[END_REF]. As initially reported for Drosophila and mouse [START_REF] Cayrou | New insights into replication origin characteristics in metazoans[END_REF] and then for human [START_REF] Langley | Genome-wide identication and characterisation of human DNA replication origins by initiation site sequencing (iniseq)[END_REF][START_REF] Besnard | Unraveling cell typespecic and reprogrammable human replication origin signatures associated with g-quadruplex consensus motifs[END_REF], G-rich sequence motifs at replication origins would have the potential to form G-quadruplexes 8 contributing to origin specication. The genome wide studies also demonstrated that active origin are associated to transcription start sites (TSS) and GC rich regions [START_REF] Langley | Genome-wide identication and characterisation of human DNA replication origins by initiation site sequencing (iniseq)[END_REF][START_REF] Besnard | Unraveling cell typespecic and reprogrammable human replication origin signatures associated with g-quadruplex consensus motifs[END_REF][START_REF] Mesner | Bubble-seq analysis of the human genome reveals distinct chromatin-mediated mechanisms for regulating early-and late-ring origins[END_REF]. In summary, if knowledge about replication origin properties made signicant progresses, the current methods cannot precisely quantify origin eciency and the relation between replication origin usage and specication of the replication timing program remains to be claried.

Analysis of replication fork directionality 1.3.1 Replication Fork directionality

The stochastic nature of replication origin locations and ring times results in the orientation of the fork movement at a given locus x to depend on the considered S phase. The Replication Fork Directionality (RFD) is dened as the dierence between the proportion of rightward (p (+) (x)) and leftward (p (-) (x)) moving forks at a given locus x: RF D(x) = p (+) (x)p (-) (x).

(1.1) So, at a strong replication origin, where replication by a leftward moving fork shifts to replication by a rightward moving fork, the RFD will changes sharply from -1 to 1, resulting in a upward jumps (Figure 1.7). If the exact position of replication initiation changes from one cell cycle to the next within a so-called initiation zone (IZ), RFD will present an upward linear slope across the IZ with an amplitude that reects the eciency of the IZ (Figure 1.7).

Correspondingly, at xed termination, the RFD jumps sharply downward and dispersed termination sites within a termination zone (TZ) lead to a linear decrease of the RFD across the TZ (Figure 1.7). Many methods were used to compute the RFD proles based on the replication fork asymmetry (Section 1.1.2). In the following section, we will detail two of them: the nucleotide compositional skew, which is an indirect method based on sequence composition, and the OK-seq method, which is a direct experimental method.

Nucleotide compositional skew

It was shown that if the replication strands are subjected to the same mutation/repair patterns, at the equilibrium, the equalities of [A]=[T] and [C]=[G] ([.] denotes the frequency of nucleotide in the DNA sequence) should be veried on each strand [START_REF] Lobry | Properties of a general model of DNA evolution under no-strand-bias conditions[END_REF]. Inversely, when the complementary strands are subject to dierent mutation patterns, the nucleotide composition of strands can become asymmetrical, which is revealed by a departure from the equilibrium compositional equalities. One of the method to quantify the strand asymmetry is to calculate directly from the sequence the compositional TA skew and GC skew:

S TA = n T -n A n T + n A , S GC = n G -n C n G + n C , (1.2)
where n A , n C , n G and n T are respectively the numbers of A, C, G and T in the analyzed sequence window.

Most bacterial genomes exhibit a GC skew, and some a TA skew [START_REF] Rocha | Universal replication biases in bacteria[END_REF][START_REF] Rocha | The replication-related organization of bacterial genomes[END_REF]. The skew proles partition the circular bacterial chromosomes in two halves, one with positive skew values and the other with negative skew values separated by two sharp jumps: an upward jump at the replication origin and a downward jump at the terminus [START_REF] Lobry | Asymmetric substitution patterns in the two DNA strands of bacteria[END_REF]. Hence, the sign of the skews correspond to the ±1 RFD proles of these chromosomes. This suggests that replication could be responsible for compositional strand asymmetry likely originating from replication coupled asymmetric mutation patterns. These skew proles became an ecient way to determine or to conrm the position of bacterial replication origins and termination sites.

The relationship observed between the compositional asymmetry and the replication pro-gram in bacteria can be generalized to eukaryotic genomes. In the human genome, as the TA and GC skews correlate [8789], the total skew S dened as the sum of the TA and GC skews:

S = S TA + S GC , (1.3)
was analyzed. It presents N shaped domains of size ranging from 0.2 to several megabases (Figure 1.6A) [START_REF] Huvet | Human gene organization driven by the coordination of replication and transcription[END_REF][START_REF] Baker | Wavelet-based method to disentangle transcription-and replicationassociated strand asymmetries in mammalian genomes[END_REF]. Based on the analogy of the bacterial case the N-domains borders, skew upward jumps were proposed to be putative replication origins [8891]. A linear decrease of the S proles was observed between these putative replication origins (Figure 1.6A). In human the S proles do not present sharp downward transitions but gradual downward slopes. Hence, the N-shaped skew prole cannot be trivially extend from the replicon model in bacteria. Because the inter-origin distance measured using DNA combing is about ∼ 50-100 kb [START_REF] Guilbaud | Evidence for sequential and increasing activation of replication origins along replication timing gradients in the human genome[END_REF], that is much smaller than the size of many detected N domains (1-3 Mb), it was suggested that some initiation and termination events must occur inside the N domains following a non-trivial random distribution which could result in the observed linear decrease of the skew. 678 N-domains bordered by 1060 putative replication origins covering ∼34% of the genome were identied in human [START_REF] Huvet | Human gene organization driven by the coordination of replication and transcription[END_REF][START_REF] Baker | Wavelet-based method to disentangle transcription-and replicationassociated strand asymmetries in mammalian genomes[END_REF].

The GC and TA skews reect the direction of the fork of replication because the leading and lagging strands experience dierent rates of nucleotide substitution [START_REF] Hyrien | Peaks cloaked in the mist: the landscape of mammalian replication origins[END_REF]. Theoretical analyses [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF] corroborated that the skew proles, in both human and bacteria, result from replication fork directionality proles by showing that under the assumption of strand-asymmetric substitution rates coupled to fork orientation, the total skew is expected to be proportional to the average fork directionality:

S(x) ∝ RF D(x).

(1.4)

The linear decrease of S proles along N-domains from positive (5' end) to negative (3' end) values would thus reect a linear decrease of the RFD with a change of sign in the middle of the N-domains i.e. replication skew domains would be the signature of large-scale gradients of replication fork polarity [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF].

Coupling between Skew and Mean Replication Timing U domains

In this section, linking replication fork polarity to MRT allows to bring another evidence of the organisation of the genome in replication domains. In Baker et al. [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF], it was shown that the replication fork polarity is related to MRT, under the central hypotheses that the replication fork velocity v is constant and that replication is bidirectional from each replication origin.

In this scenario, the replication timing prole t(x) in one cell cycle is completely specied by the position x i and the ring time Fork coming from Oi meets the fork coming from Oi+1 at termination site Ti with spacetime coordinates (yi, ui) given in Eq. [START_REF] Parker | Mechanisms and regulation of DNA replication initiation in eukaryotes[END_REF].

t i = t(x i ) of
Note that one can deduce the fork polarity in B (resp. origin and termination site locations in C) by simply taking successive derivatives of the timing profile in A (Eqs. [START_REF] Cvetic | Eukaryotic origins of DNA replication: could you please be more specic?[END_REF] and ( 20)). The fundamental hypothesis is that the replication fork velocity v is constant. Reprinted from Audit et al. 33 .

Finally, using the Dirac function to represent origin locations (x x i ) and termination sites (x y i ) (Fig. 12(c)), we obtained the following fundamental relationships:

v d dx r(x) = p(x), (19)
v d 2 dx 2 r(x) = d dx p(x) = 2 X i (x x i ) X i (x y i ) ! . (20)
In other words, we can extract, up to a multiplicative constant, the fork polarity p(x) (Fig. 12(b))

and the location of origins and termination sites (Fig. 12(c)) by simply taking successive derivative (Eqs. [START_REF] Cvetic | Eukaryotic origins of DNA replication: could you please be more specic?[END_REF] and ( 20)) of the timing profile r(x) (Fig. 12(a)). Finally, we need to rewrite the above results for application to experimental replication data that actually operate on a large number of cells (millions), from which only population statistics can be derived and have a finite spatial resolution of tens of kb or more 26,83,84,93,99,141 . Since taking the spatial derivative commutes with statistical and spatial average, we obtained:

d dx hr(x)i Cells, x = 1 v hp(x)i Cells, x , (21)
d 2 dx 2 hr(x)i Cells, x = 2 v N Ori Cells, x (x) N Ter Cells, x (x) , (22)
where h . i Cells, x stands for the average over many cells and over the spatial resolution x and

N Ori Cells, x (x) (resp. N Ter Cells, x (x)
) is the number of origins (resp. termination sites) per unit length averaged over many cells and the spatial coordinate. Note that when replication fork speed is inhomogeneous with fork speed fluctuations that do not depend on the replication timing nor on the spatial coordinate then Eqs. ( 21) and ( 22) remain valid with v standing for the average replication speed. Note that one can deduce the fork polarity in B (resp. origin and termination site locations in C) by simply taking successive derivatives of the timing prole in A. The fundamental hypothesis is that the replication fork velocity v is constant. (Adapted from Audit et al. [START_REF] Audit | Multiscale analysis of genomewide replication timing proles using a wavelet-based signal-processing algorithm[END_REF]).

The replication fork polarity prole p(x) in that cell cycle is directly deduced from the position of the origins O i and termination sites T i (Figure 1.8B), whose spatial distributions can be represented by Dirac functions of weight +1 and -1, respectively (Figure 1.8C). One can deduce the fork polarity in Figure 1.8B (and the resp. origin and termination site locations in Figure 1.8C) by simply taking successive derivatives of the timing prole in Figure 1.8A. These results can be rewritten for application to experimental replication data obtained from a large number of cells (millions), from which only population statistics can be derived with a nite spatial resolution of tens of kb or more [START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF][START_REF] Hansen | Sequencing newly replicated DNA reveals widespread plasticity in human replication timing[END_REF]9295]. Since taking the spatial derivative commutes with statistical and spatial average, one gets:

d dx MRT(x) = 1 v RFD(x) , (1.5)
d dx RFD(x) = 2(N Ori (x) -N T er (x)), (1.6)
where MRT(x), RFD(x), N Ori (x) and N T er (x) stands for the average over many cell cylces and over some spatial resolution ∆x of the replication timing, replication fork polarity and numbers of origins and termination sites per unit length.

Equations (1.4) and (1.5) reveal that the RFD provides an immediate connection between the skew S and the mean replication timing:

S(x) ∝ d dx MRT(x). (1.7)
Since the mutations responsible for the observed strand compositional asymmetries only accumulate in germline cells, this relationship was veried using, as a substitute for germline MRT, the MRT in BG02 embryonic stem cell [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF]. The skew S was indeed correlated with dMRT/dx, in the BG02 embryonic stem cells (Figure 1.6D). The signicant correlations observed in intergenic (R = 0.40, P < 10 -16 ), genic (+) (R = 0.34, P < 10 -16 ) and genic (-)

(R = 0.33, P < 10 -16 ) regions are representative of the correlations observed in an other 6 cell lines [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF]. These correlations are as important as those obtained between the dMRT/dx proles in dierent cell lines [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF], as well as those previously reported between the replication timing data themselves [START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF][START_REF] Hansen | Sequencing newly replicated DNA reveals widespread plasticity in human replication timing[END_REF][START_REF] Yae | Comparative analysis of DNA replication timing reveals conserved large-scale chromosomal architecture[END_REF]. From equation (1.7) one also deduces that the integration of the skew S is expected to generate a prole rather similar to the MRT prole. In segments where S behaves linearly, its integral is thus expected to show a parabolic prole. The integrated S function when estimated by the cumulative skew Σ ( A novel method was developed to measure directly the orientation of the fork at each locus

x of the human genome ( [3]). This method is based on the purication and sequencing of the Okazaki fragments (the method will be detailed in Chapter 3). The Okazaki fragments (Figure 1. (1.8)

The N-shaped RFD proles observed along skew N-domains and replication timing U-domains were reproduced in the RFD proles obtained by the Ok-seq method (Figure 1.10b,c,d) [3].

The remarkable similarity of RFD gradients derived from skew [START_REF] Chen | Replication-associated mutational asymmetry in the human genome[END_REF] or replication timing [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF] or determined by Ok-seq provides strong mutual validation of these data. In addition, the Introduction higher resolution of Ok-seq reveals the areas of initiation and termination of replication with much greater accuracy than with skew or MRT proles.

All together, there exist a strong correlation between the RFD and MRT proles. Under the hypothesis that the fork velocity is constant, the steep slopes at the U-domain boundaries corresponds to a high level of |RFD| whereas the at slope in the center of the U-domains corresponds to a null RFD (Figure 1.10a). In addition, the U-domains frequently coincided with N-shaped domains detected independently from the asymmetric strand composition [51 53, 90] (Figure 1.10c), proposed to reect an N-shaped pattern of RFD [START_REF] Huvet | Human gene organization driven by the coordination of replication and transcription[END_REF][START_REF] Chen | Replication-associated mutational asymmetry in the human genome[END_REF] in germline.

OK-seq fully veried the predicted RFD gradients deduced from replication timing U-domains and skew N-domains (Figure 1.10bcd). It was proposed that these linear gradients of RFD would result from cascades of replication origin ring propagating from the domain boundaries to the center [START_REF] Guilbaud | Evidence for sequential and increasing activation of replication origins along replication timing gradients in the human genome[END_REF][START_REF] Audit | Megabase replication domains along the human genome: relation to chromatin structure and genome organisation[END_REF][START_REF] Hyrien | From simple bacterial and archaeal replicons to replication N/U-domains[END_REF].

Interplay between DNA replication program and other genomic features

The modication of the replication program between cell lines can occur at regions in any stage the replication program. For example, the change of eciency of a replication origin could be associated with a change of ring time, of transcription level or of chromatin state at this locus. Hence it is of general interest to understand the link between replication initiation, replication timing, chromatin state and transcription.

Coupling between nuclear architecture and the replication program

It has been observed that within the nucleus, DNA replication is spatially organized into multiple small subnuclear domains, called replication foci [START_REF] Berezney | Heterogeneity of eukaryotic replicons, replicon clusters, and replication foci[END_REF]. The distribution pattern of replication foci reproducibly depends on S phase progression (Figure 1.5b). At the beginning of the S phase, replication foci are distributed everywhere in the nucleus. In mid of S phase, these foci tend to relocate at the nuclear membrane or at the periphery of the nucleoli. Finally, at the end of S phase, larger foci corresponding to a clustering of foci are observed at the same peripheral zones (nucleus and nucleolus). Late-stage S foci appear to be associated to heterochromatin zones (Figure 1.5b) [START_REF] Nakayasu | Mapping replicational sites in the eucaryotic cell nucleus[END_REF]. The structure of these foci and their biological function are not fully established. It has been proposed that these foci shelter several consecutive replicons, whose synchronous initiation depends on a single initiator complex. More recently, high resolution uorescence microscopy studies suggested that these foci correspond to a clustering of multiple replication complexes, which would be activated at the same time [START_REF] Baddeley | Measurement of replication structures at the nanometer scale using super-resolution light microscopy[END_REF][START_REF] Cseresnyes | Analysis of replication factories in human cells by super-resolution light microscopy[END_REF]. This was also suggested by the model proposed by Berezney et al. in 2000 [39] and taken up by Cayrou et al. in 2011 [79], according to which the replicons are clustered and synchronously replicated.

Coupling between replication and transcription

For decades, researchers have been eager to nd the connections between two fundamental genomic mechanisms, the DNA replication and the transcription of genes. Analyses in human [START_REF] Woodne | Replication timing of the human genome[END_REF] and Drosophila [START_REF] Schübeler | Genome-wide DNA replication prole for drosophila melanogaster: a link between transcription and replication timing[END_REF][START_REF] Macalpine | Coordination of replication and transcription along a drosophila chromosome[END_REF] have shown that there is a relationship between early replication timing and active gene transcription of close-by euchromatin DNA strands. However, in Drosophila, transcriptional activity has to be considered over large domains (> 100 kb) rather than at the level of individual genes to understand the impact of transcription on replication timing [START_REF] Macalpine | Coordination of replication and transcription along a drosophila chromosome[END_REF]. More recent studies have conrmed this now classic pattern of early/late replication associated with gene expression/repression but suggested that in fact other factors such as chromatin state are more likely the causally linked factors to replication timing [START_REF] Hansen | Sequencing newly replicated DNA reveals widespread plasticity in human replication timing[END_REF].

Using the SNS-seq, origin of replications were indentied in HeLa cells and it was found that most of them overlap with transcriptional regulatory elements [START_REF] Cadoret | Genome-wide studies highlight indirect links between human replication origins and gene regulation[END_REF]. Recently using the OK-seq method, Petryk et al. described 9836 (HeLa) and 5684 (GM06990) replication initiation zones corresponding to ascending segments of the RFD proles and observed that > 46% were anked by at least one active gene [3]. They also reported 9440 (HeLa) and 5715 (GM06990) large (8-2600 kb) descending segments (DS) of preferential replication termination covering 49.5% and 39.4% of the genome, individually, with 5000 (35% of the genome) shared between the two cell lines [3]. In early S phase, descending segments were associated to transcribed gene bodies and to regions of co-oriented replication and transcription. In late S phase, descending segments predominantly comprised large non-expressed DNA regions.

Association between DNA replication and epigenetic

The complex formed by DNA and its packaging proteins is named chromatin. In chromatin, ∼ 150 bp DNA segments wrap around a bead-like structure formed by an octamer of histone proteins, to form the nucleosome. Histones are the most prevalent proteins in chromatin.

There are ve types of histones: H1, H2A, H2B H3 and H4, . H3 and H4 associate together to form a dimer; two dimers of H3-H4 associate to form a tetramer which is, in turn, surrounded by two H2A and H2B dimers. The H1 protein is a linker protein (10 to 90 bp long) that completes the nucleosome. Histones carry epigenetic marks that convey a functional information by two mechanisms. First, one of the canonical histones can be replaced by a histone variant 9 .

For instance, transcriptionally active regions are enriched in H3.3 and H2AZ variants of H3 and H2A, respectively. Second, histones are formed of a globular part that constitues the nucleosome core and a exible tail that reaches outside toward the nuclear enviromnent. These tails can carry diverse covalent modications that have a functional meaning. There is a specic to regulatory proteins. The epigenetic term groups histones modications described above as well as the modications of DNA as for example cytosine methylation (Reviewed in [102104]).

Chromatin organization and its dynamics participate in essentially all DNA-templated processes including transcription and replication. Replication origin positioning is cell line specic: even if there were consensus sequence motifs specifying replication origins, an additional regulation mechanism would be needed [START_REF] Gilbert | In search of the holy replicator[END_REF]. Experiments demonstrated that origin ring depends on the chromatin environment [START_REF] Macalpine | Chromatin and DNA replication[END_REF]. Therefore, mechanisms that position and control the time of ring of origins are likely associated with epigenetics marks and linked to chromatin structure [START_REF] Méchali | Eukaryotic DNA replication origins: many choices for appropriate answers[END_REF][START_REF] Mcnairn | Epigenomic replication: linking epigenetics to DNA replication[END_REF][START_REF] Courbet | Replication fork movement sets chromatin loop size and origin choice in mammalian cells[END_REF][START_REF] Smith | Chromatin structure and replication origins: determinants of chromosome replication and nuclear organization[END_REF]. It was demonstrated that epigenetic marks and chromatin accessibility can be used to faithfully estimate the DNA replication timing proles in a cell type specic manner [START_REF] Gindin | A chromatin structure-based model accurately predicts DNA replication timing in human cells[END_REF]. From a reverse perspective, it was argued that replication timing proles can be interpreted as an epigenetic signature of a cellular dierentiation state [START_REF] Hiratani | Replication timing as an epigenetic mark[END_REF]. In both views, there is strong a coupling between the DNA replication and the epigenetic marks proles. A recent study described human chromatin as constituted of four chromatin states when analysis eleven genome wide epigenetic proles at 100 kb resolution of MRT in 5 somatic cell lines and one embryonic stem cell (ESC) [START_REF] Julienne | Human genome replication proceeds through four chromatin states[END_REF][START_REF] Julienne | Embryonic stem cell specic master replication origins at the heart of the loss of pluripotency[END_REF] For more details about the chromatin structure and its inuence on DNA replication, we refer the reader to the review in Ref. [START_REF] Macalpine | Chromatin and DNA replication[END_REF]. of the RFD proles among the 12 cell lines. We show that our classication remains stable when considering "probes" as small as 50 Mb suggesting that replication program changes are widespread over the whole genome. We identify that the regions rich in GC are more conserved among the 12 cell lines in terms of RFD proles. This suggests that there is a spatial susceptibility to replication program changes between cell lines. From these observations, we dene the regions of the genome with an unstable replication program. In chapter 4 we extend these cell line comparisons to replication timing and transcription. On the one hand, we nd the cell lines classication are similar to those obtained using the RFD proles. On the other hand, we conrm that the large late replicating domains correspond to unstable regions in terms of replication. We also address in a direct manner the association between replicative changes and transcription changes and we recover the existence of a general coupling between replicative changes and transcription changes. Comparing the expression levels between GM06990

and Raji cells, we identify 2075 genes that change expression level and their are coupled with replicative change relatively. In chapter 5 we focus on a cellular system modeling the progression of Chronic Myeloid Leukemia (CML). First, we demonstrate that by using the RFD proles we are able to follow gradual replicative changes occurring during CML progression.

Then, performing manual annotation of the eciency changes of replication initiation zones among the cell lines of the CML system, we observe that the targeted initiation zones are more frequently weakened than enhanced and that weakened initiation zones are associated with transcription repression. Chapter 6 presents a general conclusion and perspectives.

CHAPTER 2

Materials and Methods

Cell lines

We mainly report analyses concerning the characterization of replication and transcription in 12 cancer and non cancer human cell lines (Table 2.1). The 12 cell lines belongs to three types: Four lymphoid cell lines (GM06990, Raji, BL79, IARC385), four myeloid cell lines which form a Chronic Myeloid Leukemia `CML' model (TF1_GFP, TF1_BCRABL_1M, TF1_BCRABL_6M, K562) and four adherent cell lines (IMR90, HeLa, TLSE19, IB118). Table 2.1 recapitulates the main cell lines properties and cancerous status. We refer the reader to [3,[START_REF] Wu | Developmental and cancer-associated plasticity of DNA replication preferentially targets gc-poor, lowly expressed and late-replicating regions[END_REF] for the details of cell sources and culture conditions.

Adherent cell lines

IMR90 is a normal primary broblast cell line, obtained from ATCC. IMR90 cells are untransformed and presenescent cells [START_REF] Yu | Induced pluripotent stem cell lines derived from human somatic cells[END_REF]. The three other adherent cell lines are cancerous. HeLa cell line is a cancer cell line established by George Gey and commonly used in cell biology and medical research. HeLa cell line comes from a metastasis sample taken from Henrietta Lacks, an African American patient with cervical cancer who died in 1951 [START_REF] Lucey | Henrietta lacks, HeLa cells, and cell culture contamination[END_REF]. HeLa was the rst immortal cell line of human origin. Note that an immortal cell line of animal origin was created 11 years before by Wilton Earle [START_REF] Mcgehee | Johns hopkinsthe birthplace of tissue culture: the story of Ross G[END_REF]. These cells have the distinctive property of dividing indenitely, which means that they can proliferate in laboratory culture conditions.

TLSE19 and IB118 are Leiomyosarcomas (LMS) cell lines. Sarcoma refers to the malignant growth that originates from transformed cells out of mesenchymal tissues.

Lymphoid cell lines

We used two Burkitt lymphoma cell lines. The rst brief description of Burkitt Lymphoma (BL) cells was established in 1897 by doctor Albert Cook [START_REF] Hutt | Historical introduction, burkitt's lymphoma, nasopharyngeal carcinoma and kaposi's sarcoma[END_REF]. In 1958, Dennis Burkitt described in details BL cell lines [START_REF] Burkitt | Burkitt's lymphoma outside the known endemic areas of africa and new guinea[END_REF][START_REF] Magrath | African burkitt's lymphoma. history, biology, clinical features, and treatment[END_REF]. This type of Lymphoid cancer is aggressive and widespread among the population and can aect children and adults. One reason BL is of interest here is that the overexpression of Myc gene and Epstein-Barr virus (EBV) infection distinguish BL from other Lymphomas. The overexperssion of Myc in human broblast induces hyper-activation of origin ring in early S phase [START_REF] Dominguez-Sola | Non-transcriptional control of DNA replication by c-myc[END_REF][START_REF] Dominguez-Sola | Myc and the control of DNA replication[END_REF]. The First BL cell line is Raji. This cell line was established by R.J.V Pulvertaft in 1965 from a Burkitt lymphoma extract from an 11 years old male from Africa (Nigeria) [START_REF] Pulvertaft | A study of malignant tumours in nigeria by short-term tissue culture[END_REF]. Raji is the rst continuous human cell line of hematopoietic origin. The second BL cell line is BL79 duplicated from the International Agency for Research on Cancer (IARC).

We also used two control normal B lymphoblastoid cell lines (LCL): GM06690 and IARC385.

A common method to generate LCL is the transformation of peripheral B lymphocyte cells from blood by density centrifugation. While T cells are removed from the cell population, the transformation of B lymphocytes is accomplished by EBV infection. IARC385 was duplicated from IARC and was obtained from the same patient as BL79.

Myeloid cell lines

Myeloid cell lines comprise a cellular model for establishment and progression of chronic myeloid leukemia (CML), a malignant disease characterized by the Philadelphia (PH1) chromosome and the formation of the BCR_ABL1 fusion gene, whose expression is necessary and sucient for CML formation [START_REF] Daley | Induction of chronic myelogenous leukemia in mice by the p210bcr/abl gene of the philadelphia chromosome[END_REF]. TF1 is a BCR_ABL negative cell line established from a patient with erythroleukemia. TF1_GFP and TF1 BCR_ABL were obtained by transduction of green uorescent protein (GFP) or a BCR_ABL_GFP fusion gene using a murine stem cell virus (MSCV)-based retroviral vector [START_REF] Laperrousaz | Primitive CML cell expansion relies on abnormal levels of BMPs provided by the niche and bmprib overexpression[END_REF]. TF1_BCR_ABL cells were analyzed after culturing for 1 month (TF1_BCRABL_1M) or 6 months (TF1_BCRABL_6M) following transduction. K562 (the last myeloid cell line) was obtained by culturing leukemic cells of pleural eusion of a patient with PH1-positive CML in blast crisis, now called K562 [START_REF] Lozzio | Human chronic myelogenous leukemia cell-line with positive philadelphia chromosome[END_REF].

Hence, K562 is a late CML model.

Various Genomic data 2.2.1 GC Content

GC-content uctuations recapitulate the non uniform organization of gene size [START_REF] Duret | Statistical analysis of vertebrate sequences reveals that long genes are scarce in GC-rich isochores[END_REF], and gene density [START_REF] Costantini | The evolution of isochore patterns in vertebrate genomes[END_REF]. There is also a correlation between GC content and other features of the genome [START_REF] Bernardi | Misunderstandings about isochores. part 1[END_REF]. For example, GC-rich regions are poor (resp. rich) in LINE1 1 (resp. ALU 2 ) repeat se- quences. Interestingly, there is an important correlation between GC content and replication timing [START_REF] Woodne | Replication timing of the human genome[END_REF][START_REF] Costantini | Replication timing, chromosomal bands, and isochores[END_REF], early replication being correlated with high GC content.

The GC-content was computed as :

GC% = n G + n C n A + n T + n G + n C * 100 (2.1)
where n A , n T , n G and n C are the numbers of A, T, G and C nucleotides counted along the window of interest. We dened GC-content categories following the 5 isochores classication of the human genome [START_REF] Bernardi | Misunderstandings about isochores. part 1[END_REF] in light isochores L1 (GC<37%) and L2 (37%≤GC<41%) and heavy isochores H1 (41%≤GC<46 %), H2 (46%≤GC<53% ) and H3 (GC≥53%). The genome coverage of L1, L2, H1, H2, H3 was 26.5%, 32%, 24.0%, 13.1%, and 4.4%, respectively, after classication based on GC content in non-overlapping 10 kb windows.

Comparative Genomic Hybridization data

Total genomic DNA was extracted for GM06990, Raji, IARC385 or BL79 cells to determine the copy number of each loci. We represented in Figure 2.1 for each cell line the log 2 ratio of the copy number relative to a normal karyotype along of each chromosome. The red line represents the loss in copy number as we see in IARC385 cell line for the chromosomes X and 18. The blue line represents the gains in copy number as we see in chromosome 7 in IARC385 cell line. The black lines correspond to normal copy number without any gains or loss. GM06690 appears as a normal cell line where the copy number is free from amplication or deletion. 1 The long interspersed nuclear elements are part of the repeated sequences dispersed within the DNA. Their length are usually a few thousands base pairs and their exact nature varies according to the species. 2 An Alu sequence is a short DNA fragment belonging to the family of small interspersed nuclear elements or SINE. The Alu sequences are the most abundant transposable element of the human genome [START_REF] Consortium | Initial sequencing and analysis of the human genome[END_REF], distributed throughout all chromosomes, representing about 10% of the human genome [START_REF] Zhang | Regulatory roles of Alu transcript on gene expression[END_REF].

MRT computed for 7 cell lines

MRT was proled by Repli-seq [START_REF] Chen | Impact of replication timing on non-CpG and CpG substitution rates in mammalian genomes[END_REF][START_REF] Thurman | Identication of higher-order functional domains in the human ENCODE regions[END_REF], which consists in sequencing of newly replicated DNA of cells sorted by DNA content into consecutive compartments of S-phase. For GM06990, GM12878, K562, HeLaS3 and IMR90, alignment les of Repli-seq libraries (BAM les) for 6 S-phase fractions were obtained from the ENCODE project [START_REF] Hansen | Sequencing newly replicated DNA reveals widespread plasticity in human replication timing[END_REF] at http://hgdownload. cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeUwRepliSeq/. The number of aligned reads ranged from 8.5 x 10 6 in K562 to 62.8 x 10 6 in IMR90 and all reads were used in down- stream analyses. Repli-seq data for TF1_GFP and TF1_BCRABL_1M were obtained using a modied protocol based on the incorporation of EdU [START_REF] Wu | Developmental and cancer-associated plasticity of DNA replication preferentially targets gc-poor, lowly expressed and late-replicating regions[END_REF] instead of BrdU as in the original protocol [START_REF] Hansen | Sequencing newly replicated DNA reveals widespread plasticity in human replication timing[END_REF][START_REF] Chen | Impact of replication timing on non-CpG and CpG substitution rates in mammalian genomes[END_REF]. Puried nascent DNA was sequenced on illumina HiSeq platform and aligned to human genome (hg19) at the IB2C sequencing platform. Only uniquely mapping reads (186.0 x 10 6 for TF1_GFP and 143.3 x 10 6 for TF1_BCRABL_1M ) were used in downstream analyses. To compute MRT proles, tag densities in 100 kb windows were normalised and denoised as in [START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF][START_REF] Audit | Multiscale analysis of genomewide replication timing proles using a wavelet-based signal-processing algorithm[END_REF] using modied thresholds for the two EdU-based datasets (TF1_GFP and TF1_BCRABL_1M). At each locus, the distribution of replication times was obtained by normalizing the denoised tag densities of each S phase fractions by their sum.

MRT values were estimated as the mean of these distributions.

2.2.4 Quantication of RNA-seq for the 12 cell lines using TopHat For Raji, GM06990, BL79, or IARC385, total RNA was extracted from 1-10 x 10 Gene transcriptional levels for the 12 cell lines were estimated using the same computational pipeline based on the TopHat suite of softwares [START_REF] Trapnell | Dierential gene and transcript expression analysis of RNA-seq experiments with tophat and cuinks[END_REF]. TopHat (version 2.1.1) and bowtie2 (version 2.2.9) were used to align RNA-seq reads to the human genome (hg19). RNA abundance were computed using Cuinks (version 2.2.2). We fed the reference transcript annotation provided by Illumina iGenomes (ftp://igenome:G3nom3s4u$@$ussd-ftp.illumina.com/Homo_ sapiens/UCSC/hg19/Homo$_$sapiens$_$UCSC$_$hg19.tar.gz) to TopHat (-G option) and cuinks (-G option), and only considered RNA abundance estimates for genes present in the reference annotation. For each cell line we obtained, for the same set of 24,371 genes of size > 300 bp, the estimated mRNA level expressed in FPKM (Fragments Per Kilobase of exon model per Million mapped fragments) and a 95% condence interval on the FPKM ([FPKM lo ,

FPKM hi ]). We considered a gene to be expressed when FPKM ≥ 1 and we ltered out cases where FPKM hi /FPKM lo > 2 (at most 118 genes were ltered out in TF1_BCRABL_6M).

The proportion of expressed genes ranged from 43.5% in Raji to 51.3% in TF1_BCRABL_1M.

Transcription level of any region R of length l R was estimated based on the FPKM values of all the expressed genes g that overlap with the regions as:

FPKM(R) = g l g∩R l R FPKM(g), (2.2)
where l g∩R is the overlap length between g and R.

2.3 Statistical Methods

Pearson and Spearman Correlation coecients

Correlation coecients are statistical measures constructed to capture the strength of the relationship between the relative uctuations of two random variables. Their values range between -1.0 and 1.0. A correlation of +1.0 indicates a perfect positive correlation (in phase), while a correlation of -1.0 indicates a perfect negative correlation (in opposite phase). A correlation of 0.0 indicates no relationship between the uctuation of the two variables.

The most common correlation coecient is the product-moment correlation of Pearson [START_REF] Saporta | Probabilités, analyse des données et statistique[END_REF]. It is simply dened as the covariance Cov(X,Y) of two variables (X,Y) normalized by their standard deviations (σ X and σ Y ) to correct for the dependence of the covariance on the amplitude of the uctuations, thus assuring Corr(X,Y) ∈ [-1,1]:

Corr(X, Y ) = Cov(X, Y ) σ X * σ Y . (2.3)
Given a sample of n joint observations (x i , y i ) of X and Y, Corr(X,Y) is estimated as follows:

ρ XY = n i=1 (x i -x) * (y i -ȳ) n i=1 (x i -x) 2 * n i=1 (y i -ȳ) 2 , ( 2.4)
where x and ȳ are the empirical means.(x = 1 n n i=1 x i , same for ȳ).

The Spearman rank correlation coecient is dened as the Pearson correlation coecient of the ranks of X and Y [START_REF] Izenman | Modern multivariate statistical techniques, Regression, classication and manifold learning[END_REF]. Let r i be the rank of x i and s i the rank of y i , to compute the Spearman rank correlation coecient we used Eq. 2.4 replacing x i by r i and y i by s i . Pearson correlation coecient captures linear relationships between X and Y uctuations. Spearman's is more robust as it captures any monotonic relationship. If there are no repetitive data values, a perfect Spearman rank correlation occurs with +1 or -1 when the two variables are a monotonic function of the other. Note that, in this thesis, the results of the correlation analyses were qualitatively unchanged if calculated using Spearman rank-correlation coecients instead of Pearson correlations. Hence, we chose to report only results obtained using Pearson correlations.

Hierarchical ascending classication

Hierarchical ascending classication (HAC) is a simple iterative classication method [START_REF] Bruynooghe | Classication ascendante hiérarchique des grands ensembles de données: un algorithme rapide fondé sur la construction des voisinages réductibles[END_REF].

Given a partition of a set of N objects into k classes, the core idea is to join the two "closest" classes into one class and thus obtain a new partition in k-1 classes.

Starting from N classes of one object, HAC amounts to repeat this agglomeration scheme until obtaining one class of N objects. This process requires a dissimilarity measure d(o, e) between pairs of objects and an aggregation rule to derive dissimilarity between object classes. These successive groupings produce a binary classication tree that can be represented as a dendrogram, the root of which corresponds to the class grouping all the individuals. This dendrogram represents a hierarchy of partitions. We can then choose a specic partition by truncating the tree to a given level, depending on some constraints as the number of classes.

"Minimum jump" strategy or "single linkage" aggregation rule amounts to dene the dissimilarity between two classes C 1 and C 2 as their "nearest neighbors" distance:

∆(C 1 , C 2 ) = min o∈C 1 ,e∈C 2 d(o, e) (2.5)
It produces "chains" of classes, linked to each other like links in a chain, by pair of objects that are close to each other [START_REF] Sibson | SLINK: an optimally ecient algorithm for the single-link cluster method[END_REF].

In opposite, "maximum jump" or "complete linkage" aggregation rule amounts to dene the dissimilarity between to classes C 1 and C 2 as the dissimilarity between their most distant objects [START_REF] Saporta | Probabilités, analyse des données et statistique[END_REF][START_REF] Defays | An ecient algorithm for a complete link method[END_REF]

∆(C 1 , C 2 ) = max o∈C 1 ,e∈C 2 d(o, e) (2.6)
HAC algorithm does not take into account the nature of the dissimilarity function provided; it does not matter whether it is a true distance or another derived measure that may be more meaningful to the context under study. We can therefore choose a type of dissimilarity adapted to the context and the nature of the data. The most direct method to calculate dissimilarity between objects in a multidimensional Euclidean space consists in calculating the Euclidean distances. In this thesis we use the correlation distance to compute the dissimilarity between the transcription or the DNA replication proles of a pair of cell lines. It is computed as 1 minus the correlation coecient value between the pair of cell lines proles:

D(CL 1 , CL 2 ) = 1 -Corr(CL 1 , CL 2 ) (2.7)

Correlation matrices representation and comparison

For the 12 cell lines cited previously, we performed comparative analyses of their DNA replication and transcriptional programs. We illustrate the general procedure we followed using the replication fork directionality proles discussed in Chapter we computed the Pearson correlation coecients (C orr(RF D i , RF D j )) (Eq. (2.4)) between the RFD proles of all 66 pairs of cell lines (i, j). To visualize the result, we grouped the coecients in a symmetric matrix (C RF D (i, j) = Corr(RF D i , RF D j )). As can be seen on Figure 2.2A, when arbitrarily ordering the matrix lines and columns according to cell line names alphabetical order, the correlation matrix structure cannot be easily observed.

Single linkage clustering (section 2.3.2) using pairwise correlation distance (Eq. (2.7)) as the dissimilarity measure results in a hierarchical binary classication of the cell lines.

Ordering the correlation matrix lines and columns according to a cell line order respecting this classication (the lines of the dendrogram representation do not cross), we obtain a correlation matrix representation where the correlation structure between cell lines is apparent (Figure 2.2B). In the example of Figure 2.2, the 3 group classication of our 12 cell lines highlighted by this dendrogram representation, nicely matches the cell type classication described in section 2.1.

In the following chapters we use this representation strategy to illustrate the results of correlation analyses between genomic proles such as RFD, MRT (C M RT ) and gene expression level (FPKM) (C RN A-seq ) among the 12 cell lines.

We quantify the similarity between the correlation structures of two correlation matrices in two ways. First, the classication trees can be compared, in particular we ask whether these 3 group classications are identical or not. Second, we compute the Pearson correlation coecient between the vectors of 66 (non-trivial) pairwise genomic prole correlation coecients (contained in the upper triangular part excluding the diagonal of the correlation matrices). It will be referred to as the correlation coecient between two correlation matrices.

D C 2 C 1 | = |RFD(C 1 )-RFD(C 2 )|) between two cell lines (C 1 ,C 2 )
, we only consider windows with > 2000 OK-seq reads in both cell lines to guarantee that the standard deviation, using a Poisson approximation, remains < 0.023 for both RFD estimates. -G option provides to TopHat with a set of annotations of known gene models and / or transcripts.

Bioinformatic softwares

TopHat rst extracts the transcribed sequences and uses Bowtie to align the read on this virtual transcriptome. The read that maps the transcriptome is be converted to genomic mappings and merged with the new mappings and junctions of the nal output of TopHat.

no-coverage-search option disables the coverage based search for junctions, it is used to speed up computations.

The gene dierential analysis is performed using the DESeq2 [START_REF] Love | Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2[END_REF] package and its R version (3.5.1). We use Featurecount procedure to create the genes read count table across all replicates from the TopHat alignments. Dierential gene expression analysis is then performed using DESeq2 procedure on the count table.

We use Samtools (1.9) to read, view, create index and edit SAM and BAM format les [START_REF] Li | The sequence alignment/map format and SAMtools[END_REF]. To visualize the mapped read along the genome we use the Integrative Genome Viewer (IGV (2.4.8)).

Introduction

Several techniques have been used to understand the DNA replication program and describe the nature and position of replication origins in mammalian genomes [START_REF] Fachinetti | Replication termination at eukaryotic chromosomes is mediated by top2 and occurs at genomic loci containing pausing elements[END_REF][START_REF] Raghuraman | Replication dynamics of the yeast genome[END_REF]143145]. The microarray hybridization and next generation sequencing based methods allows to detect and measure the abundance of replication intermediates genome wide. These technologies include SNS-seq and Bubble-seq. Using SNS-seq method, between 50000 and 250000 potential initiation sites were identied on the human genome [START_REF] Besnard | Unraveling cell typespecic and reprogrammable human replication origin signatures associated with g-quadruplex consensus motifs[END_REF][START_REF] Cayrou | The chromatin environment shapes DNA replication origin organization and denes origin classes[END_REF][START_REF] Mukhopadhyay | Allele-specic genome-wide proling in human primary erythroblasts reveal replication program organization[END_REF]. In contrast, using Bubble-seq more than 100000 sites were identied in the human genome [START_REF] Mesner | Bubble-seq analysis of the human genome reveals distinct chromatin-mediated mechanisms for regulating early-and late-ring origins[END_REF]. The issue using these methods is that bubble trap origins align weakly with those obtained with SNS-seq method [START_REF] Langley | Genome-wide identication and characterisation of human DNA replication origins by initiation site sequencing (iniseq)[END_REF][START_REF] Cadoret | Genome-wide studies highlight indirect links between human replication origins and gene regulation[END_REF][START_REF] Karnani | Genomic study of replication initiation in human chromosomes reveals the inuence of transcription regulation and chromatin structure on origin selection[END_REF][START_REF] Mesner | Bubble-seq analysis of the human genome reveals distinct chromatin-mediated mechanisms for regulating early-and late-ring origins[END_REF] and that dierent sets of origins are detcted with these two methods.

Moreover, they do not provide the eciency of the replication origin. Interestingly, Okazaki

Fragments sequencing (Ok-seq) method, that maps initiation events along the whole-genome takes advantage of the deep-sequencing of strand-specic Okazaki fragments, to solve these issues. DNA replication is inherently asymmetric. About 50% of the genome is discontinuously replicated as Okazaki fragments, discovered in 1968 by Reji and Tsuneko Okazaki The Ok-seq method was rst established in yeast. In this rst implementation [START_REF] Mcguee | Quantitative, genome-wide analysis of eukaryotic replication initiation and termination[END_REF][START_REF] Smith | Intrinsic coupling of lagging-strand synthesis to chromatin assembly[END_REF],

the isolation of Okazaki fragments was only possible thanks to the DNA ligase I repression in a degron-tagged yeast construct. For the rst time, Smith and Whitehouse puried Okazaki fragment and showed that on average the Watson:Crick log 2 ratio of Okazaki fragment strand changed sign when crossing Active origins in yeast [START_REF] Smith | Intrinsic coupling of lagging-strand synthesis to chromatin assembly[END_REF].

In the same manner, using the deep sequencing of Okazaki fragments, the directionality of the replication fork in yeast was conrmed and modeled [START_REF] Mcguee | Quantitative, genome-wide analysis of eukaryotic replication initiation and termination[END_REF]. Very recently, Okazaki fragments sequencing (Ok-seq) was also accomplished in 2 human cells: GM06990

1 and HeLa 2 without the need for ligase inactivation [3] (for more details see section 3.2.1). This work quantied the replication fork directionality and inferred the location of replication initiation and termination zones, genome wide, at high resolution.

Here we have used novel Ok-seq 1 GM06990 is an EBV-immortalized lymphoblastoid cell line (LCL) with a near-normal karyotype 2 HeLa is an epithelial cell line from a cervix adenocarcinoma and sucient for CML formation [START_REF] Daley | Induction of chronic myelogenous leukemia in mice by the p210bcr/abl gene of the philadelphia chromosome[END_REF]. The aim of this chapter is to do a comparative analysis of the DNA replication program of these 12 human cell lines established using Ok-seq.

3.2 Proling the DNA replication program by sequencing of Okazaki fragment

Purication of Okazaki fragment

The Ok-seq method in yeast requires to use conditional lethal mutants that massively increase the Okazaki fragment abundance [START_REF] Mcguee | Quantitative, genome-wide analysis of eukaryotic replication initiation and termination[END_REF][START_REF] Smith | Intrinsic coupling of lagging-strand synthesis to chromatin assembly[END_REF]. In human a new methodology was used to isolate 3 See Chapter 5 for the denition

3).

The technique is directly applicable to mammalian cells, which readily incorporate EdU. However, wild type yeast cells are impermeable to EdU. Therefore, the application of EdU based Ok-seq to yeast requires the use of strain modied to incorporate EdU.

Estimating Replication fork Directionality using OF sequencing

Sequenced reads were identied and demultiplexed using the standard Illumina software suite and adaptor sequences were removed by Cutadapt (version 1.2.1 to 1.12). Reads >10 nt were aligned to the human reference genome (hg19) using the BWA (version 0.7.4) software with default parameters. This resulted in alignment les in 'BAM' format that were downloaded from the IMAGIF sequencing platform 4 (Figure 3.2). We only considered uniquely mapped reads and we counted identical alignments (same site and strand) as one to remove PCR duplicate reads. We computed the Replication Fork Directionality (RFD) for a given window following Eq. (3.1) as:

RF D = (R -F )/(R + F ), (3.1)
where F represents the number of Forward strand matching reads (leftward moving forks) forks and R represents the number of Reverse strand matching reads (rightward moving forks) forks within the window; RFD takes value between -1 and 1 (Figure 3.3). In the rst two top panels of Figure 3.3, we represented the number of forward reads (anti-sense forks) and 4 Centre de Recherche de Gif-sur-yvette, CNRS reverse reads (sense forks) at 10 kb scale in a 20 Mb window. We can see that the reverse and forward reads are complementary; for example, when we have an abundance of reverse reads we have a deciency of forward reads and vice-versa. This abundance of reads switching between forward and reverse determines the RFD proles (Figure 3.3).

Since the total amount of replication on both strands should be constant across the genome, normalization of the dierence between the two strands counts by the total read count (R+F) accounts for variation in read depth due to copy number or sequence bias (Figure 3.3).

As originally described in Petryk et al. [3], RFD proles display an alternation of quasilinear ascending, descending and at segments (AS, DS, FS) of varying size and slope (Figure

Application to 12 human cell lines

We used Ok-seq to prole RFD genome-wide in multiple cell lines as previously described for HeLa and GM06990 [START_REF] Ryba | Abnormal developmental control of replication timing domains in pediatric acute lymphoblastic leukemia[END_REF]. As described in the previous chapter we have three types of cell lines; Four lymphoid (GM06990, Raji, BL79, IARC385), four myeloid (TF1_GFP, TF1_BCRABL_1M, TF1_BCRABL_6M, K562) and four adherent (IMR90, HeLa, TLSE19, IB118) cell lines. From 2 (IB118, GM06990) to 6 (BL79) biological replicates (BR) were sequenced per cell line (Table 3 which increased resemblance to K562, a late CML. CML progression system will be specically analyzed in chapter 5. For the adherent cell lines, including the 2 LMS cell lines (TLSE19 and IB118), the strongest resemblance was observed between TLSE19 and IMR90. TLSE19 was only slightly more correlated to IB118 than to HeLa. IB118 was only slightly more correlated to TLSE19 than to IMR90, but was more distant to HeLa. This suggests that the cancer cell lines are classied based on their tissue of origin. Thus cell type of similar morphologies, and similar dierentiation pathways, tend to present similar RFD proles.

We repeated this genome-wide analysis for three other resolutions (100 kb, 200 kb and 1 Mb, Figure 3.7 BCD). The same classication of cell lines was obtained as for 10 Kb resolution. This showed the robustness of RFD to preserve the classication in high and low resolution. In the majority of the following RFD analyses, we chose the 10 kb resolution.

The changes in RFD proles are widespread among the genome

In this section, we are interested to know if the information is localized in specic regions or is distributed in the whole genome. Firstly, to investigate whether RFD dierences between cell lines were concentrated in particular chromosomes, we repeated the analysis shown in The minimal value of correlation is between chromosomes 16 and 22 (Figure 3.8). This result suggests that the information obtained from RFD proles is not due to specic regions, but is distributed among all chromosomes. We stratied the genome using the GC content classes described by Bernardi [START_REF] Bernardi | Isochores and the evolutionary genomics of vertebrates[END_REF]. Why did we choose the GC content? Simply because it is common to all cell lines, and is associated with MRT and transcription, so that it is a convenient substitute for genome organization

chr1 chr2 chr3 chr4 chr5 chr6 chr7 chr8 chr9 chr10 chr11 chr12 chr13 chr14 chr15 chr16 chr17 chr18 chr19 chr20 chr21 chr22 All Chromosome chr1 chr2 chr3 chr4 chr5 chr6 chr7 chr8 chr9 chr10 chr11 chr12 chr13 chr14 chr15 chr16 chr17 chr18 chr19 chr20 chr21 chr22
For each class of GC content we computed the RFD correlation matrix and performed the hierarchical classication illustrated by the dendrogram. The classication and the structure of the correlation matrix do not change (Figure 3.9) so we deduced that the dierences of RFD proles among the twelve cell lines are also widespread among the isochore families.

To more precisely assess how widespread RFD changes are, we generated a large number of random probes, 50 kb to 50 Mb in size, consisting each of 55000 randomly located 10 kb windows. For each probe, we computed its RFD correlation matrix and the correlation with the global genome correlation matrix. We then determined the statistical uctuation of the latter correlation and the probability of observing the same classication of cell lines as with the global genome (Figure 3.10). We found that a random probe size ≥5 Mb was sucient to reach a >90% correlation with the global genome RFD correlation matrix, and to observe with a probability >0.95, the same three group classications of cell lines. These latter result shows that cell line classication is not caused by outliers data point but represent the data in general.

F ∆ C L 1 RFD GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 G ∆ C L 2 RFD GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 H ∆ C H 1 RFD GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 I ∆ C H 2 RFD GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 J ∆ C H 3 RFD 0.2 0.1 0.0 0.1 0.2

Characterizing the heterogeneity of the DNA replication program

We have shown that 5 Mb probes are sucient to capture the DNA replication program changes between cell lines. However, we visualized that this is not the case in 5 Mb windows (maximally localized probes). The replicative program is in fact aected by characteristics of the genome, like the GC content. This reveals an heterogeneity along the genome in the behavior of the changes of RFD proles among cell lines.

In this section, we will determine these regions, and classify the genome into replication stable and variable regions.

Characterization of these regions, by linking to other genomes features, can help us to decipher the origin of this heterogeneity.

RFD proles are more conserved in high GC content regions

In the previous section, we showed that the global classication of RFD proles among the 12 cell lines is conserved when we stratied the genome using GC content levels. But, we can observe that the amplitude of the correlation coecient in the correlation matrix increase simultaneously with the increase of GC content levels (Figure 3.11A-E). In other words, the RFD proles among the 12 cell lines in the rich GC regions are more correlated than in poor GC regions. Furthermore, the myeloid cell lines are more correlated to lymphoid cell lines in L1, L2, H1 regions than to adherent cell lines (Figure 3.11A-C). This suggests that even in low GC content regions we can nd some high correlation between myeloid and lymphoid cell lines (Figure 3.11A). However, a better correlation between myeloid and adherent cell lines is found in high GC content regions. In contrast, even in GC-rich regions we do not have a correlation between lymphoid and adherent cell lines (Figure 3.11E). These observations can be explained by the tissue of origin of cell lines. The lymphoid and myeloid cell lines are the part of blood cells that are produced in the bone marrow, concomitantly to a specic inter-group correlation.

When the correlation coecient dierences between each GC-content class (Figure 3.11A-E) and the entire genome (Figure 3.7A) were computed (Figure 3.11F-J), most dierences were negative in L1, null in L2 and increasingly positive in H1 to H3. In other words, the RFD proles were less, equally, or more similar to each other in the L1, L2, or H1-3 fractions, respectively, than in the total genome. Therefore, RFD changes were more frequent in the GC-poor fractions of the genome. There were a few exceptions to the general trend of increasing RFD correlation in the L1<L2<H1<H2<H3 order (Figure 3.11). The correlations between HeLa and GM06990, and between IB118 and Raji, decreased rather than increased with GC content, in the order H2<H3<H1<L1<L2. In addition, the order was L1< L2<H3<H1=H2 for HeLa vs. IB118, H1=L1<H2=L2<H3 for HeLa vs. IMR90 and L1<H3<L2<H1<H2 for IARC385 vs. Raji. It was previously observed that in HeLa and GM06990, many IZs border active genes, and isolated genes expressed in only one cell type are anked by IZs only in that cell type [3]. Given that genes are enriched in GC-rich isochores [START_REF] Woodne | Replication timing of the human genome[END_REF], such transcription-dependent changes in IZ activity are expected to decrease the RFD prole correlation more strongly in GC-rich than in GC-poor regions. The comparison of IZs in GM06990 and HeLa also revealed that in addition to these gene-bordering IZs, many cell-type specic IZs are found away from active genes, in late-replicating [3] GC-poor regions.

Since GC-poor isochores form a much larger fraction of the genome than GC-rich isochores, the net density of RFD changes between HeLa and GM06990 is higher in GC-rich than in GC-poor regions, due to a predominant contribution of transcription-associated changes in GC-rich regions. Inversely, the net density of RFD changes between other cell lines is most often higher in GC-poor regions, suggesting a predominance of transcription-independent RFD changes in GC-poor regions.

Note that, low GC content level 'L1' corresponds to 26.5% and 'L2' corresponds to 32% of the genome an inversely high GC content level 'H3' corresponds to 4.4% of whole genome.

To verify if within each isochores group we have a homogeneity of correlation coecient value among the 12 cell line. We stratied the genome into 10 GC content decile groups (Appendix, Thus, the new stratication showed that the GC dependence of the replicative changes is also observed within isochores classes. We obtained the same results and interpretation as with the previous stratication (Figure 3.11).

To test if the more frequent RFD changes in GC poor regions are not due to technical noise, we analyzed in Figure 3.12 the mean value of RFD dierence correlation matrix between the global matrix of correlation of RFD proles and each matrix of correlation of RFD proles for the 5 isochores. If the RFD changes were due to noise dierences, correlations dierences should vanish when the scale of analysis is increased. However, the cell lines classication (Figure 3.7A) and the GC dependence of correlation dierences (Figure 3.11) were conserved or even enhanced when the scale of analysis was increased from 10 kb to 100 kb, 200 kb and 1Mb (Figures 3.7 and 3.12). Therefore, these observations were not due to a higher technical noise in GC-poor regions.

We have seen in this section that RFD prole changes between cell lines do not behave the same along the genome and depend in particular on the concentration of the GC. Therefore, can we nd some regions with a high GC content that have a high level of modications in RFD proles among the 12 cell lines? To address this issue we have to classify the regions according to the stability of their replication program.

Spatial heterogeneity of DNA replication program variability in 5 Mb windows

To answer the question that we addressed at the end of the previous section, we decided to divide the genome into 5 Mb windows.

Identication of regions with variable DNA replication program

To identify the regions where the RFD proles of the 12 cell lines are similar (as in Figure 3.13) or very dissimilar (as in Figure 3.14), we computed the correlation dierence matrix

Z w = 1 n Σ n i=1 x i 1 n Σ n i=1 (x i -x) 2 , ( 3.2)
where x i is the set of n=66 pairwise correlation dierences in window w. The histogram in Blue curves corresponds to the 5% most variable regions. Black curves corresponds to the whole genome. stable (resp. variable) RFD proles are more associated to GC-rich (resp. AT-rich) where 50% of the regions have a GC above (resp. under) 0.42 (resp. 0.36), thus they predominantly belong to H1, H2 and H3 isochores (resp. L1). The stable (resp. variable) replication regions are highly (resp. lowly) expressed, 50% (resp. 10%) of them have a FPKM>1 (Figure 3.17).

The median of the MRT in the variable regions is 0.7, this suggests that the changes in RFD proles tend to be in late replicating regions. However, in the stable regions, it is equal to 0.4. This shows that 50% of these regions are in early replicating regions. The method used to quantify the transcription level by FPKM and Mean Replication timing will be detailed in the following chapter. In conclusion, replication program stability appears more heterogeneous than expected (Figure 3.15) in a manner correlated to genome organization (Figure 3.17).

Mapping of DNA replication program conservation along the human genome Mean Correlation of RFD (MCR)

We introduced a novel method to map replication program conservation along the genome.

This method is based on the RFD proles of the 12 cell lines. First, we computed the mean of the 12 RFD proles in each locus x at 1 kb resolution:

D(x) = 1 n n i=1 RFD i (x), (3.3)
where x is a window of 1 kb, with RFD i (x) is the RFD proles of the i th cell line and n=12.

Then we computed for a given window (W), the correlation coecient (C 1,i (W)) between the mean RFD prole D(W) and each RFD proles for the 12 cell lines RFD i (W). In the nal step, we calculated the mean of these correlation coecients for the 12 cell lines:

M CR(W ) = 1 n n i=1 (C 1,i (W)). (3.4)
At the resolution of 500 kb, MCR values are between 0.2 and 0.9 with an average equal to 0.52 as illustrated in the Figure 3 To illustrate the consistence between the MCR approach and previous results obtained for 5 Mb probes (Figure 3.16), we computed the MCR proles at 500 kb along two previously analyzed regions with a stable (Figure 3.13) and variable (Figure 3.14) replication program.

As expected, the MCR score for the stable (Figure 3.18 B) regions is high and the MCR score for the variable (Figure 3.18 B) regions is low. We illustrate an overview of human chromosomes using MCR levels (Figure 3.20).

When comparing the MCR proles with the GC proles along chromosomes [START_REF] Bernardi | Isochores and the evolutionary genomics of vertebrates[END_REF], we observe the concordance between them, the correlation coecient is 0.82 along the whole genome. This result suggests that the low MCR . value corresponds to the low GC content regions. For example, chromosome 4 is classied as poor in GC and at the same time using the MCR score, we can see that the DNA replication program tends to be variable in this chromosome (Figure 3.20B). These results conrm the association between the low GC content and the instability of the RFD proles among the 12 cell lines. In previous sections, we have used the GC content as a classier of the genome.

We have seen that there were more replicative changes in the low GC content regions. To test whether our MCR score ranks the regions of the genome well, we tested the level of association between the MCR and the GC-content. We computed the cdf of GC content according to the MCR level and we observed that we have a high correlation between low GC regions and low MCR regions (Figure 3.21). Figure 3.20: MCR and GC are highly correlated. We represent for the chromosome 1 (panel A) and 4 (panel B) the MCR score (bottom) and the GC content (top). We can visualize some concordance between MCR at 500 kb and GC isochore [START_REF] Bernardi | Misunderstandings about isochores. part 1[END_REF][START_REF] Bernardi | Isochores and the evolutionary genomics of vertebrates[END_REF] in chromosome 1 and 4 where we have regions rich in AT (blue) corresponding to regions of low MCR value (green). Finally, in this section 3.4, using the GC content as an indicator of genome organization, we showed that the RFD proles among the 12 cell lines are more conserved in high GC content.

Then we showed that along the genome we can identify many 5 Mb windows that are classied as very variable RFD proles regions and other 5 Mb windows that are classied as very stable RFD proles regions. To increase the resolution of detection, we computed the MCR score at 500 kb. Using MCR score we can do a very precise selection of variable regions taking into account the specicity of the region according to other features of the genome. We will adapt this score in the following analysis in order to understand the source of the heterogeneity in the DNA replication program. In the following chapter we will study the association between the DNA replication program changes and other features of the genome such as replication timing and gene expression changes. data is the presence of noise. Strictly speaking, the derivative of a noisy prole is not dened and, correspondingly, the naive derivative of a prole based on the numerical dierence between two successive samples is ill-dened and numerically unstable. However, the derivative of a smoothed version of a noisy prole is well dened [START_REF] Audit | Multiscale analysis of genomewide replication timing proles using a wavelet-based signal-processing algorithm[END_REF]. In other words, the rate of signal variation has to be estimated over a suciently large number of data points. Moving average of RFD proles computed in 10 kb windows shifted by 1 kb were computed using a smooth (dierentiable) window φ with a mid-height width of 26 kb and RFD slope was dened as the increment over 20 kb of that prole:

S(x) = d dx RF D * φ(x), (3.5)
where x is the position on the genome and * stand fot the convolution operator. All slopes were nally expressed in %RFD per kb. See Figure 3.23A for a RFD slope prole computed over a 5 Mb region in GM06990. We nally recorded the position and slope of the local As can be seen in Figure 3.23B, the RFD slopes proles are strongly correlated among the 12 cell lines. In this Figure we see that the high RFD slopes in some regions are conserved among the 12 cell lines (as at position 4.6 MB). The RFD slope prole corresponding to the HeLa cell line is the noisiest due to the fact that the RFD prole is also noisy (Figure 3.23B).

IB118 has a very low RFD amplitude that varies between ±0.6 [START_REF] Wu | Determination of DNA replication program changes between cancer and normal cells by sequencing of Okazaki fragments[END_REF], and that is why the local maximum RFD slope prole is smaller than in the other cell lines. We selected 3 thresholds to detect the IZ in each cell lines, 0.5, 1 and 2. To choose the threshold values, we compute the distribution of the local maximum of RFD slopes (MaxSlope, corresponds to IZ) and of the distributions of the opposite of the slopes at the local minimum of RFD slopes (MinSlope, corresponds to Termination Zones (TZ)) (Figure 3.24). The RFD slopes proles typically vary between -2 and 2 as shown in the Figure 3.24. We observe that the distributions of the MaxSlope are more asymmetric than the distributions of the MinSlope, and can reach larger values (up to 4). The two pdfs coincide at 1 in all cell lines except HeLa and IB118 due to the issues discussed above, indicating that all MaxSlope>1 likely correspond to IZ and not to noise. Thus the threshold of 1 will be mostly used in the following analysis to select the IZ location. Note that, to select IZ with high condence we choose the threshold 2, as we see in the Figure 3.24 at this threshold, there is almost no local slope minima of the same amplitude.

0.5 was also used has a permissive threshold.

IZ conservation and cell lines classication

As dened in the previous section, we choose three thresholds to count in each cell line the number of IZ. For the rst threshold 0.5, we detected 6414 origins in IB118 and 21883 (Table 3.2) in HeLa, and in average 9000 for the other cell lines. But, as shown in the histogram, The upper triangular part of the matrix represents the observed count of common origin with slope>1%RFD per kb between pair of cell lines in 30 kb windows. The lower triangular part of the matrix represents the expected count of common IZ computed >1%RFD per kb in each cell lines. Main diagonal is the number of detected IZ of slope>1%RFD per kb in each cell lines.

at this threshold, we would also select local slope minimum (Figure 3.24). Hence we choose the threshold 1 as we discussed previously. Using this threshold the number of detected IZ decreases in the 12 cell lines. Moreover, we obtained for GM06990 (5696) and HeLa (11641) very similar numbers of origins as previously detected using the Hidden Markov Model (HMM) method [3] (GM06990: 5684 detected origins, HeLa: 9836 detected origins).

The overlapping between the two methods to select the IZ goes to 90% in GM06990 and to 84% in HeLa. So our IZ detection with this threshold is in good agreement with previous IZ detection and will be most often used in the following analyses. To be more selective and select the strongest origins for each cell lines, we can increase the threshold to 2. Now the number of detected origin decreases to 67 in IB118 and 3099 in HeLa cell line (Table 3 IZ have a mean size of 30 kb (range 6-150 kb, [3]), so to know the number of overlapping IZ among the 12 cell lines, we computed the number of common origins per pair of cell lines at 30 kb windows using the threshold 1 and 2. Note that in this interpretation, we will exclude the results obtained in HeLa and IB118. At the threshold 1, we selected and ∼68 times more common IZ than expected, thus the common IZ with a higher slope would be more conserved among the 12 cell lines. These results give us a new argument to conrm that the RFD proles are robust classiers of cell lines based on their tissue of origin.

Finally, we increased the resolution to 10 kb to test if at high resolution there were still IZ conservations among the 12 cell lines. We repeated the same analysis as done at 30 kb windows at the threshold 1. We selected between 182 common IZ (between IMR90 and K562)

and 809 common IZ (between TF1_GFP and TF1_BcrAbl_1M) (Appendix, Figure A.9).

Also we observed on average about ∼3 times more than expected common IZ among the 12 cell lines (Appendix, Figure A.10). This results conrmed that even at a high resolution the IZ are very conserved among the 12 cell lines.

Replication stable regions are dense in ecient Initiation Zones

Finally in this chapter we study the association between the GC, the MCR score and the IZ spatial distribution. We compute the density of the IZ in 5 categories of GC content (Figure regions are corresponding to a low density of IZ with a low eciency and that replication origins with the highest eciency tend to localize within replication stable, high GC content regions of this IZ density.

Conclusion

As mentioned in the introduction, the replication timing (RT) has been used for several years as a good classier of cell line. The shape of the RT prole helps to identify large domains of constant timing [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF][START_REF] Gilbert | Evaluating genome-scale approaches to eukaryotic DNA replication[END_REF][START_REF] Hiratani | Dierentiation-induced replication-timing changes are restricted to AT-rich/long interspersed nuclear element (line)-rich isochores[END_REF]. In 2008 Hiratani et al. showed that DNA replication domains were conserved between independent mESC lines [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF][START_REF] Takahashi | Genome-wide stability of the DNA replication program in single mammalian cells[END_REF]. Ryba et al. further demonstrated the power of replication proling to distinguish closely related cell types, and argued that temporal domains of replication are spatially compartmentalized in 2 structural and functional units of three-dimensional chromosomal architecture [START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF]. Therefore, classication of cell lines could be established based on this information. But the resolution of the RT proles remained low and the domains of timing rather large (100 kb -1 Mb). In this thesis, we computed the RFD proles based on the Ok-seq data, and increased the resolution to 10 kb.

We have generated novel RFD proles and have explored several approaches to compare the replication programs of twelve cancer and non-cancer cell types. To identify the changes of RFD proles among the 12 cell lines we computed the correlation matrix at 10 kb. A global, unbiased correlation approach revealed that the RFD proles clustered in three separate groups corresponding to lymphoid, myeloid and adherent cells. Interestingly, we observed that cancer-associated changes in replication did not blur their developmental origin signature.

We did not detect any convergence of the replication programs of cancer cells from dierent developmental origins (e.g. LMS vs. BLs). In contrast, within lymphoid cells, we found evidence for BL-specic replication patterns. Overall, our global correlation analyses provide evidence for recurrent replication changes along specic tumor progression pathways.

GC content (isochores) classies the genome into 5 categories (isochores: L1, L2, H1, H2, H3). Based on this partition, we computed the RFD correlation matrix of each category. We were able to reproduce a similar global classication of cell lines. Based on this, the global correlation analyses further revealed that RFD changes are widespread through the genome but more frequent in GC-poor regions. The variations of RFD proles were distributed randomly across the genome and our global classication was conserved when we choose randomly 5 Mb probes. However, heterogeneity of correlation matrix of RFD proles was observed in 5 Mb windows. Using a Z-score constructed on RFD correlations matrices, we isolated windows where the correlation values within and between the 3 types of cell lines is very strong suggesting that, in these regions the RFD proles are the same among the 12 cell lines and that they thus form stable RFD regions. We also isolated windows were replication program both within and between groups appeared variable. The characterization of selected variable RFD proles regions conrmed that those regions are associated to low GC content, low gene expression, and late replication.

The latter results encouraged us to do a global cartography of the genome to identify the zones where the RFD proles do not change among the cell lines. We computed a new score called MCR to map the conservation level of the DNA replication program at 500 kb resolution.

MCR allowed to identify the regions that are variable/stable in terms of RFD proles among the 12 cell lines. Finally, we used the derivative of the RFD proles to delineate IZ locations as regions of maximal RFD slope. We observe that the stable replication program regions are associated to high MCR value and correspond to the regions enriched in the most ecient IZ that are conserved among the 12 cell lines. Inversely, the variable DNA replication program regions are associated to low MCR values and correspond to region of low initiation zones density. In the following chapter we will study the association between replication program and other features of the genome as transcription, and we will use the MCR score and the IZ identication to study more specically the association between replicative and transcription changes.

CHAPTER 4 Association between changes in replication fork polarity proles, Mean

Replication Timing and transcription DNA serves as a template for both replication and transcription. In order to avoid collisions between the transcription machinery and the replication forks, these two processes usually take place separately in time and space [START_REF] Aguilera | Causes of genome instability[END_REF]. The dierent copies of genes do not have the same regulatory sequences, so transcription can take place at some copies while others are replicated. This is the case of the genes encoding histone H4 [START_REF] Holmes | Coordinate control and selective expression of the full complement of replication-dependent histone H4 genes in normal and cancer cells[END_REF]. Nevertheless, despite these regulatory systems, multiple studies showed that conicts between transcription and replication leads to DNA damage and double-stranded breaks [START_REF] Kim | Transcription as a source of genome instability[END_REF][START_REF] Tuduri | Topoisomerase I suppresses genomic instability by preventing interference between replication and transcription[END_REF].

Nevertheless cell cycle average transcription proles revealed that genes tend to be replicated early during S phase [START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF]. In this chapter, we will study the global association between the replication fork directionality, the replication timing and the transcription program. We will start by comparing the previous cell line classication based on the RFD proles showed in chapter 3 to the cell line classications obtained when using MRT and transcription data.

After that we will see how the conserved regions of RFD proles detected among the 12 cell lines are associated to large domains of MRT. Finally, we will quantify the link between transcriptional changes and replication modications.

Introduction

Molecularly the rst correlation between the replication timing and the transcription was established in 1980 by analyzing a dozen of cell-type specic genes. Concomitantly, the rst half of the S phase was shown to corresponds to the Giemsa between the gene expression probability and the timing of replication (r = 0.61). Interestingly, the correlation between MRT and transcription was not observed in budding yeast [START_REF] Raghuraman | Replication dynamics of the yeast genome[END_REF], suggesting that this relationship was acquired at some point during evolution [START_REF] Gilbert | Replication timing and transcriptional control: beyond cause and eect[END_REF][START_REF] Donaldson | Shaping time: chromatin structure and the DNA replication programme[END_REF] (see review [START_REF] Hiratani | DNA replication timing enters the single-cell era[END_REF]). In Drosophila, it was demonstrated that most genes have a tendency to replicate in the rst third of S phase with a high probability of being expressed independently of their replication time during this period. A strong relationship between the MRT and transcription was found only in 25% of the genes that replicate later in S phase [START_REF] Schübeler | Genome-wide DNA replication prole for drosophila melanogaster: a link between transcription and replication timing[END_REF]. The analysis of the dierentiation of mouse embroyonic stem cells (mESCs) revealed that promoters of high and low CpG density, which generally have strong and weak promoter activity, respectively, showed distinct behaviors when switching to a late replication environment, where only CpG poor promoters showed a stronger tendency towards negative transcriptional regulation [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF].

In cancer, replication timing changes were reported also to be associated with a concerted change in gene expression [START_REF] Du | Replication timing and epigenome remodelling are associated with the nature of chromosomal rearrangements in cancer[END_REF]. It was found that genes that replicate earlier are signicantly increased in expression, whereas the genes that replicate later are signicantly repressed [START_REF] Du | Replication timing and epigenome remodelling are associated with the nature of chromosomal rearrangements in cancer[END_REF].

Using Ok-seq method to compare the replication and transcription program, it was found that DNA replication initiates preferentially at the transcription start sites 3 (TSS) [3,[START_REF] Chen | Transcription shapes DNA replication initiation and termination in human cells[END_REF].

To conclude, in human genome a global correlation between the replication and transcription program was observed but some genes did not apparently these general rules.

1 G-banding is a technique used in cytogenetics to produce a visible karyotype by staining condensed chromosomes.

2 DNaseI is a specic endonuclease that hydrolyses double-stranded (ds) DNA generating tri-and/or tetraoligonucleotides having 5'-phosphate and 3'-hydroxyl termini [START_REF] Sr | 12 deoxyribonuclease i[END_REF]. 3 TSS is the rst nucleotide of DNA that is transcribed into RNA. This suggested that IARC385 may share some replication and transcription signatures of BLs, and may not be as "normal" as GM06990. Indeed, the two LCLs IARC385 and GM06990

were the least correlated cell lines within this group. These results suggest either that IARC385 was established from abnormal, but non-BL blood cells from the same patient as BL79, which seems unlikely, or that IARC385 was destabilized during immortalization, which is more plausible. Both scenarios may explain the observed genomic instability associated with replication and transcription changes reminiscent of BLs. EBV transformation occurred in vivo for the BLs but in vitro for the LCLs.

Within the adherent cells, dierent classications were obtained by RFD and RNA-seq may help to distinguish dierent types of LMS and suggested a possible dierentiation of TLSE19 and IB118, which were derived from a buttock muscle tumor and a scalp tumor, respectively. The strong correlation of the RNA-seq proles of the two LMSs may be due to the use of dierent RNA-seq methodologies for LMSs and for other cell types. Alternatively, the shared expression of cancer-specic genes, in two LMSs of dierent cellular origins, may have resulted in strong convergence of RNA-seq but not RFD proles. It is notable that the RNA-seq and RFD data better matched each other for IB118 than for TLSE19. This suggests that replication is more dissociated from transcription in TLSE19 than in IB118.

(Figures 3.7
Lymphoid, myeloid and adherent cells formed three separate MRT clusters as in RFD (Figure 4.6A). A previous analysis reported that MRT proles of nonleukemic cells are much less variant than those from leukemias [START_REF] Ryba | Abnormal developmental control of replication timing domains in pediatric acute lymphoblastic leukemia[END_REF]. In agreement, we observed a strong correlation of the MRT proles from GM06990 and GM12878, two LCLs with a near-normal karyotype (Figure 4.6A). MRT proles did not reect the CML progression, since K562 appeared more correlated to TF1-GFP than to TF1-BCRABL-1M (Figure 4.6A). This cannot be explained by variations in MRT methodology since both TF1 derivatives were proled by the same method.

In summary, RFD proles suggested the existence of CML-specic replication changes, not necessarily reected in MRT proles. To further investigate this, we compared the cumulative distributions of MRT changes between TF1-GFP, TF1-BCRABL-1M and K562 to those observed between distinct LCLs (GM06990 vs. GM12878) and distinct cell types (K562, GM06990, IMR90) (Figure 4.

6B). The proportion of MRT changes with

|∆M RT | > 0.2 after one month of BCR-ABL1 expression in TF1 cells (1.13%) was slightly higher than between LCLs (0.3%), but much lower than between distinct cell types (K562 vs. IMR90: 22.9%; K562 vs. GM06990: 12.8%; IMR90 vs. GM06990: 21.8%). In contrast, the proportion of |∆M RT | > 0.2 changes between TF1-GFP and K562 (9.45%) was more similar to that previously observed between LCLs and leukemia. We concluded that RFD changes induced in our model for CML establishment and early progression are not accompanied by large-scale shifts in MRT, but that such shifts may appear during progression to late CML (K562).

In summary, globally similar results were obtained by RNA-seq, but divergences between RNA-seq and RFD classications were also observed. These results suggested that recurrent replication changes occur in specic tumor types but that the closeness of their connection with transcription may depend on cell type.

Cell line classication based on RFD, FPKM and MRT proles per chromosome

To investigate whether RNA-seq and MRT dierences between cell lines were concentrated in particular regions of the genome, we repeated the analyses shown in Figures 4.5 The results obtained for the entire genome were recapitulated for each chromosome, with minor exceptions detailed below. We recap that lymphoid, myeloid and adherent cells formed three separate RFD clusters, except that for chromosome 19, HeLa clustered with myeloid rather than adherent cells, and for chromosome 16, HeLa was equally distant from both groups. Distances within and between groups were conserved, except that for chromosome 22, myeloid cells were closer to adherent than to lymphoid cells, and for chromosome 17, IARC385, rather than GM06990, was the most distant among lymphoid cells. As to RNA-seq data, the clustering of HeLa with myeloid rather than adherent cells was again observed for each chromosome, although for chromosomes 10 and 11

HeLa appeared equidistant from myeloid and adherent cells. In addition, IMR90 clustered with HeLa and myeloid cells rather than adherent cells for chromosomes 6 and 20, and was equidistant from both groups for chromosomes 13, 14, 18 and 21. As to MRT data, the global classication was reproduced for each chromosome except that the branching point of HeLa and K562 was somewhat unstable. The correlation coecients were again smaller by RFD than by RNA-seq or MRT for each chromosome, but this was more pronounced for GC-poor chromosomes. vealed a more subtle deviation of replication and transcription changes. For K562 vs. HeLa, the order was L1<L2<H1<H2=H3 by RFD but L2<L1<H3<H2<H1 by RNA-seq. For K562 vs. IMR90, the order was L1<L2<H1<H2<H3 by RFD but L1<L2< H3<H2=H1 by RNA-seq.

GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 F ∆ C L 1 RNA-seq GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 G ∆ C L 2 RNA-seq GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 H ∆C H 1 RNA-seq GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 I ∆ C H 2 RNA-seq GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 J ∆ C H 3
Other inter-group comparisons also showed a dierent ordering of RNA-seq and RFD correlations, but in these cases we cannot exclude an eect of the dierent RNA-seq methodologies.

To summarize, contrarily to RFD changes the RNA-seq proles did not reveal a consistent tendency for transcription changes to increase or to decrease with GC content. This suggests that replication changes were at least partly dissociated from transcription changes, to an extent that depended on the cellular context.

Association between Replication Fork Polarity conservation and Mean Replication Timing

It was suggested that 50% of the genome contains large late and early domains of replication timing. The other half of MRT proles is made of U-shaped domains [START_REF] Audit | Megabase replication domains along the human genome: relation to chromatin structure and genome organisation[END_REF]. A very high association was observed between the nucleotide compositional Skew N-domains and the replication timing U domains. About half of detected U-domains are shared between dierent cell lines [START_REF] Julienne | Embryonic stem cell specic master replication origins at the heart of the loss of pluripotency[END_REF] and this is a hallmark of the robustness of MRT proles. In this thesis, MRT proles shown in Figure 4.2 help us to study the link between regions of stable and variable RFD proles among the 12 cell lines and the constant timing regions. We observed that the high MCR level is correlated to early replication but some surrounding regions are also associated to late replication. This conrmed the previous results showed in Figure 4.10AC. We further analyzed the RFD slopes proles among the 12 cell lines in this region (Figure 4.11C). We observed that this region is associated to a high density of conserved IZ (corresponding to slope maxima >1) among the 12 cell lines. In contrast, in the region between 270 and 274 Mb, which corresponds to a late replicating domains, we do not observed any IZ among the 12 cell lines. This is also true in large domains of late replication presented in Figure 4.12. Here we chose a region where we have two consecutive large late domains.

We observe a high concordance between the M RT and the MCR scores, large domains of low MCR value are highly linked to large late replicating domains.

These is compatible with the result shown in Figure 4.10BC. This two large domains are associated to a very low density of IZ among the 12 cell lines. Thus we concluded that the early replicating domains are likely linked to a high density of conserved IZ among the 12 cell lines and the late replicating domains are associated to low density of cell line specic IZ.

To further study the conservation of IZ among the cell lines, we look for some correlations between the replication IZ, the MRT and MCR proles in a pair of lymphoid cell lines (GM06990 and Raji). We selected all the RFD slope maxima (SM) > 1%RFD per kb in GM06990. Then we estimated the RFD slope change at these loci in Raji. In To conclude we described two categories of regions with a conserved replication program.

The rst one corresponds to the regions that replicate early, dense in conserved IZ have a stable RFD proles among the 12 cell lines and covers about ∼ 8% of the genome. The second one corresponds to regions that replicate late with low density of cell line specic IZ and variable RFD proles and that covers about ∼ 8% of the genome. Hence, it is important to study why the replication IZ eciency is changing among cell lines, and to correlate this change to other features of the genome such as gene expression changes.

Quantication of the links between the replication program modications and transcriptional changes

In all recent genome-wide studies, where the replication origins were identied, a link between transcription and replication origins position was detected. It appears that current studies revealed some origins of replication lying within R-loops 4 . These tri-catenary structures correspond to the hybridization of an RNA on one of the two strands of DNA [START_REF] Chedin | R-loop structures are novel, conserved, functional elements in mammalian genomes[END_REF][START_REF] Sanz | Prevalent, dynamic, and conserved R-loop structures associate with specic epigenomic signatures in mammals[END_REF])).

These results raised the hypothesis that the origins of replication activated at the beginning of S phase could be localized nearby of R-loops. These results echo an article published in 4 R-loop is a three-stranded nucleic acid structure, composed of a DNA, RNA hybrid and the non-template single-stranded DNA Blue dots represent the early replicating regions according to the MRT of GM06990 at 100 kb (resp. stable replicating regions according to MCR at 500 kb). Orange dots represent the late replicating regions according to the MRT of GM06990 at 100 kb (resp. variable replicating regions according to MCR at 500 kb). Grey dots represent the rest of the genome. 2015 that raised the question of the presence of replication origins at the level of R-loops [START_REF] Lombraña | R-loops and initiation of DNA replication in human cells: a missing link?[END_REF]. Moreover, the primary results obtained in two studies conducted in human cells [START_REF] Cadoret | Genome-wide studies highlight indirect links between human replication origins and gene regulation[END_REF] and on mouse [START_REF] Sequeira-Mendes | Transcription initiation activity sets replication origin eciency in mammalian cells[END_REF] showed that gene rich regions are also rich in replication origins and that a large part of the replication origins are at the level of transcription promoters. In this case, replication initiation sites could coincide with Transcription Starting Sites (TSS).

Cayrou et al. [START_REF] Cayrou | Genome-scale analysis of metazoan replication origins reveals their organization in specic but exible sites dened by conserved features[END_REF], showed that the strong representation of TSS in mouse genome is in fact a consequence of their association with CpG islands. In Drosophila, among the activated origins located within GC rich regions, the TSS were concomitant with gene transcription, transcriptional regulation by the mythylation of lysine 36 of histone H3 (H3K36me1 ) (reviewed in [START_REF] Wagner | Understanding the language of lys36 methylation at histone h3[END_REF]), and RNA polymerase II binding [START_REF] Comoglio | High-resolution proling of drosophila replication start sites reveals a DNA shape and chromatin signature of metazoan origins[END_REF].

Here, we used the MCR scores to delineate regions that are variable in terms of DNA replication program among the 12 cell lines. In Figure 4.15A, we represent the cumulative distribution functions of the average of FPKM values on the 12 cell lines according to their MCR levels. We found that regions with low MCR are associated to regions which are almost not expressed. About half of the regions with high MCR level (>0.57) are associated to expressed regions (Figure 4.15A). Moreover, the expression dierences between IARC385 and TLSE19, where we observed the most RFD changes, tend to conrm that the RFD changes are associated to the largest relative FPKM changes (Figure 4.15B). However, we observed that the highly expressed genes preferentially lay within stable RFD proles regions (Figure 4.15A).

Coupling between relative transcriptional changes and replication changes

In order to analyze mean expression changes between 2 cell lines, we can either consider the relative transcriptional dierences (4.2) or the absolute transcriptional dierences (4.1) between pairs of cell lines.

∆ absolute (F P KM cell 2 cell 1 ) =< |F P KM cell 1 -F P KM cell 2 | >, (4.1)
∆ relative (F P KM cell 2 cell 1 ) =< |log 10

F P KM cell 1 F P KM cell 2 | >, (4.2)
where < . > stands for the average over all windows. We compile all pairwise average absolute expression dierences (Eq. (4.1)) in a global absolute expression dierence matrix (Figure 4.16). Then we computed the same matrix but for each MCR levels from RFD variable regions (MCR<0.42, wich corresponds to 8% of all genes) to RFD conserved regions (MCR>0.57, which corresponds to 37.5% of all genes) (Figure 4.17A-D). Finally, we compute the dierences between the average absolute dierence matrix in each MCR level and the global absolute dierence matrix (Figure 4.17E-H). The blue matrix in Figure 4.17E reveals that we have less absolute transcriptional changes in the variable RFD regions except when associating GM06990 and the other lymphoid cell lines. In contrast, the red matrix in Figure 4.17H

illustrates that we have more absolute transcriptional changes in stable RFD regions. This suggests that we have a coupling between the absolute transcriptional changes and replication changes. Complementary, we computed the global relative expression dierence (Eq. (4.2)) matrix among the 12 cell lines then we follow the same procedure as done previously with the absolute change (Figure 4.18). The blue matrix in Figure 4.19H reveals that we have less relative transcriptional changes in the stable RFD regions. In contrast the red matrix in Figure 4.19E illustrates that we have more relative transcriptional changes in variable RFD proles.

Thus, we observe a coupling between relative transcriptional changes and replication changes but in the opposite direction to the coupling observed with the absolute dierences. Note that the results obtained using the absolute FPKM dierences are likely due to the level of gene expression, since the stable RFD prole regions are more highly expressed than the variable RFD prole regions and absolute expression changes are correlated to the average expression level. To summarize the results, we computed the average of the dierence matrices (obtained over the 66 independent terms) for both absolute and relative expression changes (Appendix, 4.20: Association between gene expression and RFD proles in GM06990 and Raji. (A) RFD proles and (B) FPKM proles of GM06990 (blue) and Raji (green) computed in 10 kb windows along a 5 Mb region of chromosome 2. (C) Cdf of the MCR values at the dierentially expressed gene promoters (red) and at all promoters (black). (D) Cdf of the length of gene expressed in both GM06990 and Raji (orange), all genes (black), dierentially expressed genes (red), expressed in Raji only (green) and expressed in GM06990 only (silver). (E) Counts of MCR level at the promoters of gene expressed in GM06990 and Raji (orange), all promoters (black), dierentially expressed genes (red), genes expressed in Raji only (green) and genes expressed in GM06990 only (silver).

GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 E ∆ D G1 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 F ∆ D G2 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 G ∆ D G3 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 H ∆ D G4
GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 E ∆ D G1 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 F ∆ D G2 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 G ∆ D G3 GM06990 Raji BL79 IARC385 TF1_GFP TF1_BcrAbl_1M TF1_BcrAbl_6M K562 Hela IMR90 TLSE19 IB118 H ∆ D G4

No systematic coupling between initiation zone eciency changes and FPKM changes

To complete the analysis, we increased the resolution of the analysis by comparing the gene expression level in pairs of cell lines. We selected two lymphoid cell lines: GM06990 and Raji. Using Feature count package (section 2.3.5), we counted mapped RNA-seq reads for each genes in each cell lines as the input to DESeq2 analysis. Then, we run the DESeq2 to detect the signicant gene expression changes between the two cell lines. The advantage of applying the analysis using the RNA-seq data from the 12 cell lines instead of the data from the 2 considered cell lines is that the pairwise comparison in this case takes into account the variability of genes in all cell lines. The result of DESeq2 showed Thus, we have a conrmation of the previous classication results based on RFD and RNA-seq.

We only considered the signicantly dierentially expressed genes having a log 2 fold change of gene expression >2. We chose this high threshold value to be specic in the selection of dierentially expressed genes. Under this assumption, we selected 2045 signicant gene expression changes between Raji and GM06990. First, we observed that we selected the longest length genes (Figure 4.20D). Moreover, we observed more genes specically expressed in GM06990 than in Raji (

Link between gene expression change and IZ at promoters

We further asked whether the local replication program changes are coupled with local gene expression changes? We thus correlated the IZ eciency changes to changes of the RFD slopes to the output of DESeq2. For each dierentially expressed genes, we compared the RFD slope values at the promoter in GM06990 (Normal) and Raji (Cancer), 2 cell lines of the same type.

The histogram in Figure 4.21A shows that the RFD slope changes follow on average the gene expression changes as expected. The blue line shows that we have a larger RFD slope in GM06990 than in Raji when the FPKM GM 06990 is greater than the FPKM Raji . This result conrms the previous nding in [3,[START_REF] Chen | Transcription shapes DNA replication initiation and termination in human cells[END_REF], where the replication initiation eciency increased near the active TSS increasing co-orientation of replication and transcription at gene 5' ends.

However, covariation between FPKM and RFD slope changes at dierentially expressed gene promoter was not systematic as negative ∆ slope Raji GM 06990 were observed for FPKM GM 06990 > FPKM Raji and positive genes for FPKM GM 06690 < FPKM Raji . We also chose to compare 2 cell lines with very close RFD proles, TF_GFP and TF1_BcrAbl_1M (Figure 4.21B).

We do not observe a clear co-variation between replication and transcription, likely because RFD changes are too few. Finally, we decided to manually annotate the association between dierentially expressed genes and RFD proles changes in GM06990 and Raji. We did not put any threshold on expression changes for this analysis, so we used 6060 instead of 2045 dierentially expressed genes. Only the rst 4 chromosome were annotated. We visually analyzed the RFD proles in GM06990 and Raji ±100 kb of each. Then we dened 4 types of observed changes. Type I when we observed a specic IZ for GM06990 within the ±100 kb window (Figure 4.22A). Type II when we observed to a specic IZ for Raji (Figure 4.22B).

Type III when we observed for the 2 cell lines an IZ that may change position and eciency or have the same position and eciency (Figure 4.22C). And the case that does not belongs to any of these 3 types are called 'Other' type. We nd that 38% (653/1722) of the signicant expression genes changes are associated to type III (Table 4.1) where an IZ is observed in the two cell lines, and 25% (427/1722) of the signicant expression genes changes are associated to type I and II, where an IZ is observed in the cell line of greater promoter activity only. The slope dierences between GM06990 and Raji (∆slope Raji GM 06990 = slopeGM06990-slopeRaji) at the promoter of signicantly expressed genes with FPKMGM06990 > FPKMRaji (green), with FPKMGM06990 < FPKMRaji (black), with no signicant gene expression change in green. (B) Same analysis between TF1_GFP and TF1_BcrAbl_1M. Table 4.1: Manual annotation of initiation zones eciency changes at promoter regions of dierentially expressed genes between GM06990 and Raji in for chromosomes 1, 2, 3 and 4. First column represents the chromosome. Second one represents genes orientation. The third represents the count of GM06990 specic IZ (Type I), the fourth represent the count of Raji specic IZ (Type II), the fth column for the common IZ detected in the two cell lines that may change eciency or shift position (Type III), nally the last column represents the remaining cases (Other).

A B FPKM GM06990 > FPKM Raji FPKM GM06990 < FPKM Raji No diff gene No diff gene FPKM GFP < FPKM 1M FPKM GFP > FPKM 1M

Link between IZ eciency changes and transcription changes

The reciprocal analysis considering the regions of changing RFD rather than those with changing FPKM promoter conrmed that we have a decoupling between replication and transcription changes in some regions of the genome (Figure 4.23). We selected the location where the initiation eciency of GM06990 are higher than in Raji. We rst computed the slope dierence prole as slope dierences : ∆slope Raji GM 06990 =slope GM 06990 -slope Raji in 10 kb non overlapping windows. We then selected the position of the largest ∆slope Raji GM 06990 >1.5%RFD per kb. Then to each selected position we associated the FPKM value of the nearest transcribed gene (FPKM>1) for both cell lines. Finally, we computed the expression dierences (∆F P KM Raji GM 06990 = F P KM GM 06990 -F P KM Raji ) and the relative expression changes (log 10 F P KM GM 06990 F P KM Raji ). We observed that the distribution of both expression dierences and relative expression changes between GM06990 and Raji is shifted to the right of the corresponding distributions over the complete gene set. We thus observed a global increase of gene expression at region of enhanced replication initiation activity, as expected. In the same manner, we inversed the analysis by taking into account that the regions where the RFD slope in Raji is higher than in GM06990. We observed that the distributions of FPKM dierences and relative expression changes between GM06990 and Raji is shifted to the left of the global distribution, corroborating our rst observation. We repeated the same analysis for the two lymphoid cell lines (TF1_GFP and TF1_BcrAbl_1M) (Figure 4.23CD), but did not observe expression modications at the selected regions of high replication initiation eciency change. This is likely due to the small expression modications between these two related cell lines.

In order to conrm these results, we then focused on replication IZ identied as local maxima of the RFD slope proles. Analysis of 452 IZ strongly ecient in lymphoblastoid cell line GM06990 (>2%RFD per kb) and signicantly less ecient in Burkitt These results conrm previously reported correlations between transcription and replication changes but underline that it is not systematic and that the source of the link between IZ eciency variations and transcriptional changes during cell dierentiation remains to be fully understood.

Conclusion

It is widely accepted that the replication of the human genome proceeds in a time-ordered manner, correlated with parameters as gene activity, chromatin structure, and nuclear location. However, the evidence supporting this view depends on the scale of analysis. For example, chromosome replication bands have shown a correlation between R bands 5 and gene density. Fine scale, genome-wide analysis revealed that these correlation patterns were not as strong. The early and late replication domains detected by Repli-seq method correlate very well with the subnuclear A (euchromatic, active) and B (heterochromatic, inactive) compartments dened by Hi-C [START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF][START_REF] Yae | Comparative analysis of DNA replication timing reveals conserved large-scale chromosomal architecture[END_REF]. This observation provided a better understanding of the complex association between replication timing and transcription. Moreover, other analyses using the replication timing of individual genes revealed a strong correlation between transcription and early replication [178182]. Expressed genes replicate at early stages of S phase, whereas not transcribed genes replicate at late stages of S phase [183185]. However, it was also found that same active genes are additionally found in the late replicating regions and that in the early replicating regions there are dormant genes [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF][START_REF] Farkash-Amar | Global organization of replication time zones of the mouse genome[END_REF]. This result was also observed recently by Gilbert team [START_REF] Dileep | Rapid irreversible transcriptional reprogramming in human stem cells accompanied by discordance between replication timing and chromatin compartment[END_REF], as they found that the average changes in transcription are coordinated with replication timing changes but can be anti-correlated for individual gene. Experiments in chicken cells illustrated this complex connection between transcription and replication timing [START_REF] Hassan-Zadeh | Usf binding sequences from the hs4 insulator element impose early replication timing on a vertebrate replicator[END_REF]. It was demonstrated that origins anked by the transcription factor binding site USF1 6 (for "Upstream Stimulatory Factor") contains HS4.

The HS4 insulator has the capacity to impose a shift to earlier replication. This adjustment in replication timing does not require that transcriptional action be related with it despite the fact that it is enhanced when a gene is transcribed close-by.

In a rst part of this chapter, we extended the global unbiased correlation approach to demonstrate that MRT proles clustered in three separate groups corresponding to lymphoid, myeloid and adherent cells. Similar results were obtained by RNA-seq, except that HeLa clustered with myeloid instead of adherent cells. Therefore, cancer-associated changes in replication do not blur their developmental origin signature, although changes in gene expression may sometimes do so. It is notable that the two LMSs are more correlated to each other by RNA-seq than by RFD, which suggested that their cell of origin may be dierent and that the selection for a tumor phenotype may have resulted in a stronger convergence of their transcription than their replication program. We did not detect any evidence for convergence of the replication or transcription programs of cancer cells from dierent developmental origins (e.g. LMS vs. BLs). In contrast, within lymphoid cells, we nd evidence for BL-specic replication and transcription patterns. Overall, our global correlation analyses provide evidence for the existence of recurrent replication and transcription changes along specic tumor progression pathways. Interestingly, RFD changes induced by BCR-ABL1 expression 5 R bands are guanine-cytosine rich, and adenine-thymine rich regions are more easily denatured by heat.

in early CML were not associated with large-scale MRT switches comparable to those observed between early and late CML or previously reported between leukemias and control

LCLs [START_REF] Ryba | Abnormal developmental control of replication timing domains in pediatric acute lymphoblastic leukemia[END_REF]. The global correlation analyses further revealed that RFD changes between cell lines are widespread through the genome but more frequent in GC-poor regions. In contrast, RNA-seq changes do not vary uniformly with GC content. These results strengthen the notion that replication changes are dissociated from transcription changes, to an extent that specically depends on the compared cell types.

In a second part, we aimed to understand and characterize the stable and variable regions of RFD proles among the 12 cell lines detected using the MCR proles. We found that the largest late replicating regions are associated to large domains of low MCR proles (i.e. variable regions of RFD). In contrast the largest early replicating regions are associated to large domains of high MCR proles (i.e. stable regions of RFD). Then, we determined replication initiation zones locations based on the RFD slope proles. Early replication and RFD stability was associated with a high density of the most ecient IZ that tend to be conserved between cell lines. In contrast, regions of late replication associated to variable RFD proles were characterized by a low density of IZ with comparatively smaller eciency that tend to be cell line specic. Furthermore, we quantied directly the association between the replicative changes and transcription changes. We found that there is a coupling between RFD change and relative transcriptional change on average. These conclusions conrmed the results in the literature [3,[START_REF] Chen | Transcription shapes DNA replication initiation and termination in human cells[END_REF]. We validated this analysis by doing pairwise comparison of gene expression and RFD slope proles. We observed that the signicant transcriptional changes are not systematically associated to a large RFD changes. In the same manner, we found that the initiation zones eciency changes are not systematically associated to large transcription modications. We concluded that there is no causal association between transcriptional changes and initiation zones eciency changes.

CHAPTER 5 Application in Chronic Myeloid Leukemia System

To characterize cancer cells, many tools have been developed. These include the determination of phenotypic properties, using cell culture [START_REF] Kim | Three-dimensional in vitro tissue culture models of breast cancera review[END_REF]. Other studies focused on the specic dierences in genes expression patterns between the normal and cancer cell lines [START_REF] Ertel | Pathway-specic dierences between tumor cell lines and normal and tumor tissue cells[END_REF]. Ross et al. [START_REF] Ross | A comparison of gene expression signatures from breast tumors and breast tissue derived cell lines[END_REF] explored the relationships between patterns of gene expression in breast tumor derived cell lines and those from clinical tumor specimens. It was also shown that these cell lines and tumor samples have distinct gene expression patterns [START_REF] Ross | A comparison of gene expression signatures from breast tumors and breast tissue derived cell lines[END_REF]. The unsupervised classication of this cell lines led to the separation of the breast tumor derived cell lines from those from clinical tumor specimens. Along the same line, sets of genes responsible for the dierences between solid tumors and cell lines in their response to cancer have been identied by Serial Analysis of Gene Expression [START_REF] Stein | A serial analysis of gene expression (SAGE) database analysis of chemosensitivity: comparing solid tumors with cell lines and comparing solid tumors from dierent tissue origins[END_REF]. Recently, a comparative analysis of the replication timing proles in six human cancer cell lines was performed and dened common replication timing regions between them [START_REF] Hadjadj | Characterization of the replication timing program of 6 human model cell lines[END_REF].

In order to study how the DNA replication program diers between cancerous and noncancerous cell lines, we will compared the RFD proles of 4 myeloid cell lines that form a system of evolution of a Chronic Myeloid Leukemia cancer tumor with time.

Introduction

The research is underway to determine the metastatic behavior of cancer cells in vivo. Cancer is linked to the uncontrolled proliferation 1 of cells resulting in disruption of tissue homeostasis 2 . Cancer is rst and foremost a DNA disease and the environment is associated with this process.

The development of a tumor occurs in successive stages, a single alteration of the DNA is not enough. The genome of cancer cells becomes more distinct from that of normal cells as cancer progresses [START_REF] Zink | Nuclear structure in cancer cells[END_REF]. Unlike normal cells, cancer cells are characterized by genetic instability, which is related to a deciency of the genome repair systems. This instability allows the accumulation of DNA alterations [194] that can be from genetic or epigenetic source. These genetic abnormalities may be due to the intervention of exogenous or endogenous factors. Alterations include alterations during replication (accidental) and alterations after replication (physical, chemical, biological, etc.) [4]. These genetic instabilities cause so called DNA Replication Stress (RS). RS is a term that extensively characterizes obstacles in DNA replication and by large incorporates slowing down and breakdown of DNA replication forks [START_REF] Zeman | Causes and consequences of replication stress[END_REF][START_REF] Berti | Replication stress: getting back on track[END_REF]. The meaning of replication stress is constantly developing and hard to dene accurately. RS emerges from a wide range of sources and has various ramications in the cell [START_REF] Zeman | Causes and consequences of replication stress[END_REF]. It can be a re-replication or a silencing of the replication origins [START_REF] Zeman | Causes and consequences of replication stress[END_REF]. In addition, there are other extrinsic factors that damage DNA. Intracellular problems also cause RS such as sites that are dicult to replicate, chromatin accessibility, or collisions between the replication and transcription machines [START_REF] Zeman | Causes and consequences of replication stress[END_REF][START_REF] Mazouzi | DNA replication stress: causes, resolution and disease[END_REF][START_REF] Tubbs | Endogenous DNA damage as a source of genomic instability in cancer[END_REF].

Oncogene expression can induce RS and trigger DNA damage from the earliest tumorigenesis stages [199203]. In precancerous lesions, RS induces a DNA damage response (DDR) that can trigger senescence or apoptosis. Tumorigenesis proceeds when the DDR is downregulated (e.g. by p53 mutation), favoring cell proliferation with genome instability [199203].

Oncogenes have been proposed to trigger RS by multiple mechanisms: reduced or increased origin ring, exhaustion of limiting nucleotides or replication factors, increased transcription and replication-transcription conicts. For example, in Xenopus egg extracts, in which no transcription takes place, addition of recombinant Myc increases origin ring, fork stalling, and DNA breakage in a manner dependent on Cdc45, a limiting origin ring factor, and these eects are recapitulated by addition of recombinant Cdc45 alone [START_REF] Dominguez-Sola | Non-transcriptional control of DNA replication by c-myc[END_REF][START_REF] Srinivasan | Cdc45 is a critical eector of myc-dependent DNA replication stress[END_REF]. In contrast, overexpression of HRASv12 in cultured cells stimulates RNA synthesis and RS in a manner dependent on the TATA-binding protein (TBP), a general transcription factor, and these effects are recapitulated by overexpression of TBP alone; increased origin ring seems to be a consequence rather than a cause of RS in this case [START_REF] Kotsantis | Increased global transcription activity as a mechanism of replication stress in cancer[END_REF]. Recently, a novel nascent DNA mapping assay was used to show that overexpression of Cyclin E1 or MYC, which shortens G1 phase, induces novel intragenic origins, normally erased by transcription during G1, that 1 Cell proliferation is the process that leads to an increase in the number of cells. It is determined by balancing scissions and cell loss by cell death or dierentiation. Increase cell proliferation in tumors. 2 Process involved in maintaining a xed internal state in specic tissues of the body, including control of cell proliferation and cell death and control of metabolic function. of the observed ber. Recent techniques as nanopore sequencing [START_REF] Hennion | Mapping DNA replication with nanopore sequencing[END_REF][START_REF] Mueller | Capturing the dynamics of genome replication on individual ultra-long nanopore sequence reads[END_REF] will soon allow the direct detection of stalled forks on genomic positions. However, in this chapter we will use the RFD proles to compare the eciency of the initiation zones [3] in a cell line system of cancer proportion.

Chronic Myeloid Leukemia System

To study the association between cancer and replication program, we choose to compare the IZ eciency change in a Chronic Myeloid Leukeamia (CML) [START_REF] Lozzio | Human chronic myelogenous leukemia cell-line with positive philadelphia chromosome[END_REF][START_REF] Kurzrock | The molecular genetics of philadelphia chromosomepositive leukemias[END_REF] We computed the dierences in RFD proles in non overlapping 200 kb windows between the TF1_GFP and TF1_BcrAbl_1M, and selected respectively the 1, 5, 10 and 20% largest RFD changes. We associated to these windows the GC content, MRT in TF1_BcrAbl_1M and gene expression level of TF1_BcrAbl_1M. We showed that the largest changes in the replication fork directionality are associated with poor GC content, late replication timing and low expression genes (Figure 5.3). Same results were obtained for the largest RFD changes between TF1_BcrAbl_1M and TF1_BcrAbl_6M. These observation are consistent with the results reported in chapters 3 and 4. The largest RFD changes between 1 month and 6 months of BCR_ABL1 expression are observed in GC-poor, lowly-expressed and late replicating regions. Cumulative distribution functions (cdf) of GC content (A), transcription in TF1_BrcAbl_1M (B) and MRT in TF1_BrcAbl_1M (C) computed in non-overlapping 200 kb windows of the 22 autosomes. Cdfs were determined for all windows (all, black) or selected windows with the 20% (blue), 10% (green), 5% (yellow) and 1% (red) largest RFD changes between TF1_BrcAbl_1M and TF1_BrcAbl_6M.

As discussed in chapter 3, RFD proles allow to classify the cell lines in accordance to their tissue of origin. Here, we zoom in the global RFD correlation matrix to clearly see the relation between the TF1 cell lines and K562 (Figure 5.4A). The correlation values between K562 and the TF1 cell lines increase in accordance with early CML progression. This result conrms the robustness of the RFD as a good identicator of the replicative change in this cancer progression system. It could be related to specic origin activation or silencing.

We reached similar conclusions the similar using RNA-seq (Figure 5.4B). In contrast, the relationship between K562 and TF1 cell lines is inverted when using MRT. Note that RFD proles are more discriminant than the RNA-seq, since the dierences between the correlation values of TF1_GFP and TF1_BCRABL_6M with K562 are about 0.2 for RFD proles and < 0.05 for RNA-seq proles (Figure 5.4). However, in the two cases the evolution of tumor cancer with time is clearly illustrated.

The interplay between the replication origin activity and the genome instability in yeast uncovered a role for the organization of the DNA replication in delimiting the genetic instability [START_REF] Gómez-Escoda | The organization of genome duplication is a critical determinant of the landscape of genome maintenance[END_REF]. Moreover a thousand of mutations were identied and found to be heterogeneous within and between tumors [START_REF] Mcgranahan | Clonal heterogeneity and tumor evolution: past, present, and the future[END_REF]. They could be generated by replication stress associated instability [START_REF] Negrini | Genomic instabilityan evolving hallmark of cancer[END_REF]. Hence, we hypothesized that the detected dierences in RFD correlation values are related to the eciency changes of the replication Initiation Zones. We thus compared the IZ among the three TF1 cell lines as a possible signature of RS. were conrmed (n=6) or enhanced (n=14) whereas the majority (57%) was silenced (n=24)

or weakened (n=3) in K562. In summary, there was a signicant overall tendency for the activity changes observed in K562 at these 773 IZs to occur in the same direction as at Step 2, reminiscent of the tendency of Step 2 changes to occur in the same direction as Step 1.

However, a stronger tendency to silencing was observed in K562 than at Step 2 (Table 5.1).

Complementarily, we decided to classify the detected IZ eciency changes according to their Mean Correlation Replication score (MCR) (section 3.4.3). We computed the cumulative distribution function of MCR at 500 kb for each innovation in each step. We observed that the weakened IZ eciency changes are more associated to high MCR value in step 2 ∼50%

of weakened IZ eciency changes between TF1_BCR_ABL_1M and TF1_BCR_ABL_6M are higher than 0. 4). These results suggest that the Enhanced and Weakened IZ eciency changes might be related to dierent replication modifying processes than New Silenced IZ events. In summary, BCR_ABL1 expression during early CML progression changed replication predominantly in GC-poor, lowly expressed and late replicating regions. The targeted IZs were more frequently weakened than enhanced. Targeted IZs in early CML tended to further change activity in the same direction at later tumour progression stages. Therefore, BCR_ABL1 had a long-lasting action on IZs but the direction of the change depended on the targeted region. These results suggest a potential mechanism for generating RS and genome instability independently of transcription by perturbed replication of GC-poor, late-replicating gene deserts.

Weakened initiation zones between 1 month and 6 months of BCR-

ABL1 expression are associated with transcription repression.

Many studies suggest a strong correlation between replication timing and transcription in drosophila [START_REF] Schübeler | Genome-wide DNA replication prole for drosophila melanogaster: a link between transcription and replication timing[END_REF] and mammals [START_REF] Woodne | Replication timing of the human genome[END_REF][START_REF] Woodne | Replication timing of human chromosome 6[END_REF]. They found a strong correlation between DNA replication at the beginning of the S phase and transcriptional activity. Replicating changes are associated with transcription changes for weak promoters more than strong promoters [START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF].

In the Figure 5 Therefore, although the largest RFD changes in early replicating regions were associated with transcription repression, changes in RNA expression did not reciprocally predict RFD changes. 6M between TF1_BCRABL_6M and TF1_BCRABL_1M for windows with FPKM (geometric) mean > 0. Results in green, red and purple correspond to the 5% of the regions plotted in blue with the largest absolute FPKM ratio,the largest FPKM ratio and the lowest FPKM ratio, respectively. Chapter V

Conclusion

IZ eciency changes initiated by BCR_ABL1 expression are accentuated during prolonged BCR_ABL1 expression in TF1 and furthermore in the late CML cell line K562. BCR_ABL1 therefore has a long-lasting action on IZ eciency, in a direction that depends on each IZ but is more often repressive. Changes in IZ eciency induced by BCR_ABL1 mostly occur in GC-poor, lowly expressed, late-replicating DNA regions and are predominantly repressive. These results suggest that BCR_ABL1 expression may generate RS independently of transcription by disturbing the replication of GC-poor, late-replicating gene deserts. This is in apparent contrast with other oncogenes such as MYC or Cyclin E which, by shortening G1 phase, induce early ring of intragenic origins normally erased by transcription during G1 phase [START_REF] Macheret | Intragenic origins due to short g1 phases underlie oncogene-induced DNA replication stress[END_REF]. However, late-replicating, gene-poor regions of the genome were not interrogated in the latter study.

The BCR_ABL1 kinase aects a wide range of intracellular signaling pathways [START_REF] Rowley | Chromosomal translocations and genome rearrangements in cancer[END_REF] which might directly or indirectly modulate origin ring. The origin licensing factor Cdc6 is upregulated in a BCR_ABL1-dependent manner in primary CML and K562 cells [START_REF] Zhang | Deregulated expression of Cdc6 as BCR/ABL-dependent survival factor in chronic myeloid leukemia cells[END_REF], which may explain the increased activity of some origins. The DNA damage response (DDR) pathway is activated in chronic CML [START_REF] Takagi | Atm-dependent DNA damage-response pathway as a determinant in chronic myelogenous leukemia[END_REF], and DDR activation is known to repress late-ring origins.

Even though IZ changes were predominantly repressive, both repression and activation events were observed, and the direction of most changes was conserved during CML progression, suggesting that dierent classes of late IZs have opposite responses to BCR_ABL1. It is possible that the down regulation of some late IZs by BCR_ABL1 indirectly stimulates other late IZs due to increased availability of limiting origin ring factors. It has previously been argued that replication of the human genome involves a superposition of ecient initiation at "master" IZs detected in RFD proles followed by more random, cryptic initiation between them [3,[START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF][START_REF] Guilbaud | Evidence for sequential and increasing activation of replication origins along replication timing gradients in the human genome[END_REF][START_REF] Hyrien | From simple bacterial and archaeal replicons to replication N/U-domains[END_REF]. While we only detect silencing of master IZs, it is possible that BCR_ABL1 has a much broader eect on dispersed initiation.

CHAPTER 6 Conclusion and perspectives

We have analyzed novel RFD, gene expression and MRT datasets and have explored several approaches to compare the replication and transcription programs of twelve cancer and non-cancer cell types. A global, unbiased correlation approach revealed that the RFD (Figure Therefore, cancer-associated changes in replication do not blur their developmental origin signature, although changes in gene expression may sometimes do so. It is notable that the two LMSs are more correlated to each other by RNA-seq than by RFD, which suggests that their cell of origin may be dierent and that the selection for a tumour phenotype may have resulted in a stronger convergence of their transcription than their replication program. We did not detect any evidence for convergence of the replication or transcription programs of cancer cells from dierent developmental origins (e.g. LMS vs. BLs). In contrast, within lymphoid cells, we found evidence for BL-specic replication and transcription patterns.

Furthermore, within myeloid cells, we found that expression of the BCR_ABL1 oncogene in TF1 cells, which models the establishment and early progression of CML, altered the RFD and transcription proles of TF1 cells in a manner that increased their resemblance to K562, a late CML cell line (Figure 5.4). Overall, our global correlation analyses provide evidence for the existence of recurrent replication and transcription changes along specic tumour progression pathways. Interestingly, the RFD changes induced by BCR_ABL1 expression in early CML were not associated with large-scale MRT switches comparable to those observed between early and late CML or previously reported between leukemias and control LCLs [START_REF] Takagi | Atm-dependent DNA damage-response pathway as a determinant in chronic myelogenous leukemia[END_REF]. The global correlation analyses further revealed that RFD changes between cell lines are widespread through the genome but more frequent in GC-poor regions (Figure 3 Previous works have demonstrated a covariation between the relative transcription changes and the replication program changes from a gene promoter perspective. Transcriptional changes were rst identied and, RFD proles [3,[START_REF] Chen | Transcription shapes DNA replication initiation and termination in human cells[END_REF] or MRT proles [START_REF] Dileep | Rapid irreversible transcriptional reprogramming in human stem cells accompanied by discordance between replication timing and chromatin compartment[END_REF] changes were then questioned at these loci. Similarly, we selected dierentially expressed genes between GM06990 and Raji and analyzed the change of replication initiation potential at their promoters by comparing their RFD slopes. Higher initiation potential was on average associated to higher expression level but not systematically (Figure 4.21). We conrmed this analysis by manually annotating the presence/absence of an IZ ±100 kb of the promoter of all dierentially expressed genes between GM06990 and Raji (Table 4.1). We found that, for 38% of dierentially expressed genes, an IZ was present in proximity to their promoter in the two cell lines and that, for 25% of them, an IZ was observed in the cell line of greater promoter activity only, as expected for a direct link between promoter activity and replication initiation. These results suggested that the presence of a replicative change at promoters of changing activity is in fact limited. This conrmed previous studies where it was shown that the coupling between the replication and transcription changes is not systematic for all genes.

Originally, we studied the reciprocal association at replication IZ (IZ point of view). First, we observed that the regions with a large initiation potential increase between the GM06990 and Raji are associated to increased transcription (Figure 4.23). Focusing specically on detected IZ, we selected 452 highly ecient GM06990 IZs presenting a signicant weakening in Raji. We observed that not all the IZ eciency changes are associated to a large FPKM change, as (i) 119 IZ (26.3%) are not associated to transcription in either cell lines (closest expressed protein coding gene at least 50kb away from the IZ), (ii) 46 (10.2%) were associated to gene activation in Raji and (iii) 171 (37.8%) were associated to not dierential expressed genes. Thus, we detected a global correlation between transcriptional changes and replication program changes as in the other studies [3,[START_REF] Chen | Transcription shapes DNA replication initiation and termination in human cells[END_REF][START_REF] Dileep | Rapid irreversible transcriptional reprogramming in human stem cells accompanied by discordance between replication timing and chromatin compartment[END_REF] and we quantied a signicant extent of decoupling when considering individual genes and individual IZ, underlying that other factors remain to be identied.

OK-seq method has will undoubtedly be generalised to other systems. In the hypothesis that Okazaki fragment (OF) has a constant life time all along the genome, it is expected that Ok-seq coverage is homogeneous (up to mappability or other biases). We observed that it is in fact heterogeneous, with coverage peaks sometimes corresponding to high eciency IZ Recently high-resolution Hi-C interaction maps [START_REF] Lieberman-Aiden | Comprehensive mapping of long-range interactions reveals folding principles of the human genome[END_REF] have revealed that the human genome is organized into distinct units, the so-called Topologically Associating Domains (TADs), where genomic interactions are strong within a domain and depleted between domains [START_REF] Guelen | Domain organization of human chromosomes revealed by mapping of nuclear lamina interactions[END_REF]. TADs were suggested to be a stable property of the human genome as they appeared to be conserved between dierent cell lines [START_REF] Guelen | Domain organization of human chromosomes revealed by mapping of nuclear lamina interactions[END_REF]. The comparative analysis of replication timing data and Hi-C correlation matrix in human revealed that early and late replicating loci occur in separated compartments of open and closed chromatin, respectively [START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF][START_REF] Lieberman-Aiden | Comprehensive mapping of long-range interactions reveals folding principles of the human genome[END_REF]. Moreover, the boundaries of the cell type specic replication timing domains often coincide with the boundaries of insulated compartments of chromatin interaction [START_REF] Pope | Topologically associating domains are stable units of replication-timing regulation[END_REF] and recent study showed that the TAD borders are enriched in initiations zones [3]. These results question whether the cell line specic initiation zones observed in the late replicating regions are associated with cell line specic TAD borders. q q q q q q q q q q q q q q q q q q q q q q q q q q q q q -80 -40 PC1: 42% variance PC2: 33% variance group q q q q q q q q q q q q q q q q q q q q q q q q q q q q q BL79:R1 First column represent the short name of cell lines. second column represent the state of expressed gene, third column represent the procedure of culture. "M": Myc was over-expression; "m": Myc was not over-expressed; "E": EBNA was expressed; "e": EBNA was not expressed; "l": long-term (≥ 1 month) culture; "s": short-term (1 day) culture.

12 1

 12 and copying mechanism . . . . . . . . . . . . .

Figure 1 . 1 :

 11 Figure 1.1: Semiconservative model of DNA replication.DNA replication, also known as DNA duplication, is the process in which DNA is synthesized in vivo. This mechanism makes it possible to obtain, from one DNA molecule, two molecules identical to the initial molecule, for distribution to the two daughter cells during cell division. Because nucleotide A will successfully pair only with T, and G with C, each strand of a DNA double helix labeled here as the S strand and its complementary S' strand can serve as a template to specify the sequence of nucleotides in its complementary strand. In this way, both strands of a DNA double helix can be copied precisely. In this model, the two new molecules are formed by one original (template) strand and a newly synthesized strand, this model is termed semiconservative. (Adapted from Essential Cell Biology fourth edition[START_REF] Alberts | Essential cell biology[END_REF]).

Figure 1 . 3 :

 13 Figure 1.3: Replication of eukaryotic chromosomes. The electron micrograph shows DNA replicating along a portion of chromosome of an early y embryo. The particles visible along the DNA are nucleosomes, structures made of ∼ 150 bp of DNA wrapped around a protein complexe. The drawings represent the same portion of chromosome as it might appear at dierent times during replication: (1) prior to the time of the micrograph; (2) at the time of the micrograph; (3) after the time of the micrograph and after the merger of the two rightmost replication bubbles despicted in (2). The orange lines represent the two parental DNA strands; the red lines represent the newly synthesized DNA strands. (Adapted from Essential Cell Biology fourth edition [8]).

Figure 1 . 4 :

 14 Figure 1.4: Mechanism of DNA replication at replication forks. The DNA molecules being replicated contain Y-shaped junctions, called replication forks. Two replication forks are formed at each origin of replication. At each fork, a replication machine moves along the DNA, separating both strands of the double helix and using each strand as a template to form a new daughter strand. The two forks move away from the origin in opposite directions, unfolding the double helix of the DNA and replicating the DNA as it goes. Replication of DNA in bacterial and eukaryotic chromosomes is therefore described as bidirectional. The upper diagram shows two replication forks moving in opposite directions; the lower diagram shows the same forks a short time later. Because both of the new strands at a replication fork are synthesized in the 5 to 3 direction, the lagging strand of DNA must be made initially as a series of short DNA strands, which are later joined together. To replicate the lagging strand, DNA polymerase uses a backstitching mechanism: it synthesizes short pieces of DNA (called Okazaki fragments) in the 5 to 3 direction and backward to fork movement along the template strand (toward the fork) before synthesizing the next fragment. (Adapted from Essential Cell Biology fourth edition [8]).

Figure 1 . 5 :

 15 Figure 1.5: Relationship between isochore properties and replication timing regulation, subnuclear positioning, and transcription. (a) Isochores with unusual sequence properties are subject to replication timing regulation. The mammalian genome is partitioned into isochores with dierent GC content, Long Interspersed Nuclear Element (LINE) composition, and gene density, which are generally correlated. Isochores that are high in GC and gene density but low in LINE density are replicated early in S phase, whereas the alternate extremes are replicated late.Isochores with intermediate or mixed sequence features are frequently subject to replication timing regulation during dierentiation (speculatively labeled "Facultative heterochromatin")[START_REF] Hiratani | Global reorganization of replication domains during embryonic stem cell dierentiation[END_REF][START_REF] Ryba | Evolutionarily conserved replication timing proles predict long-range chromatin interactions and distinguish closely related cell types[END_REF]. (b) Changes in replication timing that traverse the middle of S phase accompany changes in subnuclear position and transcriptional potential. Replication early in S phase (patterns I and II) takes place within the interior euchromatic compartment, whereas replication later in S phase (patterns III, IV, and V) takes place at the nuclear periphery (pattern III), at the nucleolar periphery ( pattern III), and at internal blocks of heterochromatin (patterns IV and V)[START_REF] Berezney | Heterogeneity of eukaryotic replicons, replicon clusters, and replication foci[END_REF][START_REF] Mcnairn | Epigenomic replication: linking epigenetics to DNA replication[END_REF].(Adapted from Rhind et al.[START_REF] Rhind | DNA replication timing[END_REF]).

Figure 1 . 6 :

 16 Figure 1.6: Comparing skew S = T -A T +A + G-CG+C and mean replication timing (MRT).(A) S prole (Equation (1.3)) along a 11.4 Mb long fragment of human chromosome 10 that contains 6 skew N-domains (horizontal black bars) bordered by 7 putative replication origins O1 to O7. Each dot corresponds to the skew calculated for a window of 1 kb of repeat-masked sequence. The colors correspond to intergenic (black), (+) genes (red) and (-) genes (blue). (B) Corresponding cumulative skew prole Σ obtained by cumulative addition of S-values along the sequence. (C) MRT proles from early, 0 to late, 1 for BG02 (green), K562 (red) and GM06990 (blue) cell lines. (D) Correlations between S and dMRT/dx, in BG02 (100 kb windows) along the 22 human autosomes; colors as in (A). (E) Average dMRT/dx proles (± SEM) in 663 skew N-domains after rescaling their length L to unity; colors as in (C). (Adapted from Baker et al.[START_REF] Baker | Replication fork polarity gradients revealed by megabase-sized U-shaped replication timing domains in human cell lines[END_REF]).

  Cadoret et al., using SNS-seq, identied 283 origins of replication in HeLa cells over 1% of the genome and demonstrated that the density of active origins correlated with GC content [77]. By identifying 150 new origins of replication in this same cellular model, Karnani et al. revealed the inuence of transcriptional regulation and chromatin structure on the selection of origins [78]. Moreover, the work of Karnani et al. suggested that in HeLa cells, initiation of replication takes place in AT-rich regions, near transcription initiation sites, and in open

Figure 1 . 7 :

 17 Figure 1.7: Replication fork directionality. (Top) Sketch representing replication bubble locations in 4 distinct S phases at a xed origin, a dispersed initiation zone, a xed termination site and a dispersed termination zone, from left to right respectively. The red and blue lines inside these structures indicate the Okazaki fragments synthesized during the progression of each fork. These populations of Okazaki fragments are randomly sampled from asynchronous growth cell populations in the OK-seq protocol [3]. (bottom) RFD proles expected at the corresponding loci. (Adapted from Petryk et al. [3].)

  the activated bidirectional replication origins O i , termination time and position T i being dened as the time and positions when the replication fork coming from O i meets the replication fork coming from O i+1 (Figure 1.8A).

Figure 12 .

 12 Figure 12. Modeling the spatio-temporal replication program in a single cell. (A) Replication timing r(x), (B) replication fork polarity p(x) (Eq. (2)) and (C) spatial location of replication origins (upward arrows) and termination sites (downward arrows). Oi = (xi, ti) corresponds to the origin i positioned at location xi and firing at time ti.Fork coming from Oi meets the fork coming from Oi+1 at termination site Ti with spacetime coordinates (yi, ui) given in Eq.[START_REF] Parker | Mechanisms and regulation of DNA replication initiation in eukaryotes[END_REF]. Note that one can deduce the fork polarity in B (resp. origin and termination site locations in C) by simply taking successive derivatives of the timing profile in A (Eqs.[START_REF] Cvetic | Eukaryotic origins of DNA replication: could you please be more specic?[END_REF] and (20)). The fundamental hypothesis is that the replication fork velocity v is constant. Reprinted from Audit et al.33 .

27

 27

Figure 1 . 8 :

 18 Figure 1.8: Modeling the spatio-temporal replication program in a single cell. (A) Replication timing r(x), (B) replication fork polarity p(x) and (C) spatial location of replication origins (upward arrows) and termination sites (downward arrows). Oi = (xi, ti) corresponds to the origin i positioned at location xi and ring at time ti. Fork coming from Oi meets the fork coming from Oi+1 at termination site Ti with space-time coordinates (yi, ui).Note that one can deduce the fork polarity in B (resp. origin and termination site locations in C) by simply taking successive derivatives of the timing prole in A. The fundamental hypothesis is that the replication fork velocity v is constant. (Adapted from Audit et al.[START_REF] Audit | Multiscale analysis of genomewide replication timing proles using a wavelet-based signal-processing algorithm[END_REF]).

Figure 1 . 9 :

 19 Figure 1.9: Principle of determination of RFD based on the strand of Okazaki fragment. Schematic leftward (blue) and rightward (red) forks with Okazaki fragments. RFD is computed as the dierence between the proportions of rightward-and leftward-moving forks. (Adapted from Petryk et al 2016 [3]).

Figure 1 .

 1 6B) along N-domains, indeed displays a U-shaped (parabolic) prole likely corresponding the MRT prole in the germline. Remarkably, N-domains often correspond to genome regions where the MRT in the BG02 embryonic stem cells and other cell lines is U-shaped (Figure 1.6C). Correspondingly, averaging dMRT/dx over skew N-domains in BG02 and other cell lines resulted in N-shape proles (Figure 1.6E) further validating equation (1.7) and that replication timing U-domains and skew N-domains both result from underlying N-shaped RFD proles.

Figure 1 .

 1 Figure 1.10: OK-Seq corroborates the replicative organization of N/U-domains. (a) U-domains of replication timing (top) are predicted to show an N-shaped RFD prole (bottom). The steep slopes at U-domain borders predict a high |RFD| whereas the at slope at the central, late-replicating valley predicts a null RFD [51, 52, 56]. Blue and red arrows, expected proportions of Okazaki fragments from either strand across zones of predominant initiation (IZ) at borders and predominant termination (TZ) at center of U-domain. Dark line, expected RFD prole. The mathematical relationship dMRT/dx=RFD/v implies that RFD increases (Ascendant Segment 'AS') when the timing prole is convex (d 2 MRT/dx 2 >0) and decreases (Descendant segment 'DS') when it is concave (d 2 MRT/dx 2 <0; note that the time axis is oriented from top to bottom). (b) HeLa (orange) and GM06990 (teal blue) replication timing proles showing four adjacent U-domains[51]; (c) nucleotide compositional skew prole showing skew N-domains [90] matching replication timing U-domains; (d) N-shaped RFD proles determined by OK-Seq matching N-domains. (Adapted from Petryk et al. [3]).

  4) are produced by the discontinuous replication of the lagging strand in the opposite orientation to the fork movement. So by determining the strand of origin (Watson or Crick) of Okazaki fragments, we can determine the direction of the replication fork (Figure 1.9). Counting the number W (resp. C) of Okazaki fragments on the Watson strand (resp. Crick strand), we can deduce the proportion of forks replicating in the 5' → 3' (p (+) = C/(C + W )) and the 3' → 5' direction (p (-) = W/(C + W )) and, in turn the replication fork directionality (Eq. (1.1)): RF D(x) = C -W C + W .

Introduction 1 . 5

 15 Thesis objectives Understanding how DNA replication program is adapted dynamically to normal or pathological dierentiation cases still remain the focus of great interest. This thesis work proposes an original characterization of the plasticity of the DNA replication program during development and in association with cancer progression, based on the proling of the high-resolution replication prole obtained for 12 cell lines by purication and sequencing of Okazaki fragments (Xia Wu, Hyrien Team, IBENS, Paris). From these experimental methods, we determine the replication fork directionality (RFD) genomic proles at a resolution of 10 kilobase (kb), and we performed a comparative analysis of these proles with other markers of the organization of genome activity (eg gene transcription, replication times). The manuscript is organized in six chapters. The rst one is the current introduction. Chapter 2 is the description of selected cell lines and the statistical methods (Hierarchical classication, Spearman/Pearson correlation) and informatics tools (Python packages) used in the following chapters. The results are reported in the Chapters 3, 4 and 5. In the Chapter 3, we describe in details the Ok-seq biological experiment and the RFD computing methods. Then we test the robustness of the RFD proles to classify the cell lines based on their tissue of origin by doing a comparison

Figure 2 . 1 :

 21 Figure 2.1: CGH array analysis of GM06990, Raji, BL79, and IARC385. The log2 ratio of copy number is plotted along the length of each chromosome. Red and blue horizontal segments indicate copy number losses and gains, respectively.

Figure 2 . 2 :

 22 Figure 2.2: RFD correlation matrix. (A) RFD correlation matrix where each square correspond to a correlation coecient value between two cell lines, cell lines are ordered by alphabetic of order. (B) (top) Dendrogram of the HAC, (bottom) same matrix as in A with a cell line order compatible with HAC.

2. 3 . 4

 34 Computational analyses of RFD, RNA-seq and MRT proles All analyses are restricted to the 22 autosomes to avoid artefacts due to the XX or XY karyotypes of the studied cell lines. RFD prole Pearson correlation C RF D between a pair of cell lines was computed using non-overlapping windows with > 100 OK-seq reads in both cell lines. RNA-seq correlation C RN A was computed between log 10 (FPKM) of genes expressed in both cell lines. MRT correlation C M RT is computed between non-overlapping 100 kb windows with a valid MRT estimate. Finally, when analysing GC-content, FPKM and MRT distributions in 200 kb windows depending on the change of RFD (|δRF

  Computation and Figures are implemented in python (version 2.7.13) using numpy (version 1.11.3), scipy (version 0.18.3), pandas (version 0.23.4) and matplotlib (version 2.0.0) scientic computing modules. We use Jupyter notebook as an interactive platform when needed. To compute the transcription expression level for the RNA-seq datasets we use TopHat, a fast splice junction mapper for RNA-Seq reads. It aligns RNA-Seq reads to mammalian-sized genomes using the high-throughput short read aligner Bowtie, and then analyze the mapping results to identify splice junctions between exon (Tophat 2.1.1) [135]. After alignment of the RNA-seq reads to the genome we used cuinks (2.2.2) to compute the FPKM (Fragments Per Kilobase of exon model per Million mapped fragments) values. Tophat command lines: tophat-p8-G. ./ genes. gtf--no-coverage-search-o. ./ Tophat_ out. ./ Sequence/ genomelist_ of_ Fastqfile_ per_ cell -p option indicates how many threads TopHat must use to align reads.

[ 148 ,

 148 149]. The Okazaki fragments are short and are formed on the lagging strand during DNA replication. The Okazaki fragments, as illustrated in Figure 1.4, are in the opposite direction to the replication fork progression. Okazaki fragment synthesis on the lagging strand necessitates the repeated production of single-stranded DNA and polymerization in the opposite direction to fork progression. So sequencing the Okazaki fragments yields a direct measure of the directionality of replication fork progression location of strong ring origins and termination sites (Figure 1.7). Okazaki fragments are short, between 1000 and 2000 nucleotides long in prokaryotes and are roughly 100 to 200 nucleotides long in eukaryotes [150], and they have a short life time, estimated to few seconds. Hence, there is a low amount of Okazaki fragments in wild type cells. Hence, the biggest challenges have been to nd a method to isolate the Okazaki fragments.

Figure 3 . 1 :

 31 Figure 3.1: Ok-seq protocol. Okazaki fragment purication schema (Section 3.2.1) (see Petryk et al. [3] for more information about the method).

51 HWIFigure 3 . 2 :

 5132 Figure 3.2: Example of the alignment le content for one of TLSE19 biological replicates From this le we extracted and counted the Forward strand mapping reads (ag= 0) and Reverse strand mapping reads (ag=16) to compute the RFD proles. Flag (blue square), mapping chromosome and position are in column 2, 3 and 4 respectively. When using bwa to perform read alignment, optional eld XT:A:U indicates uniquely mapping reads.

Figure 3 . 3 :

 33 Figure 3.3: Genome wide proling of replication fork directionality. From Top to bottom : number of Okazaki fragments reads mapping the reverse (R) and Forward (F) strand, corresponding RFD prole (Eq. 3.1) and total read number (R+F), computed within 10 kb non overlapping windows along 20 Mb fragment of chromosome 3 for GM06990.

Figure 3 . 4 :

 34 Figure 3.4: Extraction information from Replication Fork directionality proles. RFD prole along a 6 Mb fragment of chromosome 5 where some zones of preferential initiation (IZ, orange), preferential termination (TZ, green) and unidirectional fork progression (thick red (resp. blue) arrows: sense (resp. antisense) fork progression) are marked; the pair of red and blue thin dashed arrows mark a region of null RFD where forks progress in equal proportion in the two directions.

3. 4 )

 4 . RFD often reached values >0.9 or <0.9, indicating nearly complete purity of Okazaki fragments. A positive (negative) RFD value indicates that forks move predominantly rightward (leftward) in the cell population (Figure 3.4). Ascending Segment (AS) and Descending Segment (DS) therefore represent zones of predominant replication initiation (IZ) and termination (TZ), respectively (Figure 1.7). The amplitude of the RFD shift across each zone reects its net initiation or termination eciency. FS of high RFD are unidirectionally replicating regions. FS of null RFD, sometimes found in the middle of a TZ, are replicated equally often in both directions, presumably by random initiation and termination. A few exemplary IZs, TZs and FSs in GM06990 are illustrated for a 6 Mb segment of chromosome 5 on Figure 3.4.

Figure 3 . 5 :

 35 Figure 3.5: Replicated RFD proling are highly coherent. Classication of RFD proles obtained in each biological replicate experiment (Table2.1). Classication of cell lines using RFD proles at 50 Kb (Section 2.3.3). Pearson correlation coecient values are color-coded from blue (0.4) to red (1) using the color bar on the rights. Each master-square in the diagonal correspond to the correlation value between biological replicates of the same cell line.

Figure 3 . 6 :Figure 3 . 7 :

 3637 Figure 3.6: Genome-wide proling of replication fork directionality. (A) Replication fork directionality (RFD) along a 20 Mb segment of chromosome 3 in 4 lymphoid (top), 4 myeloid (middle) and 4 adherent (bottom) cell lines, as indicated in the top right corners of each panel.

Figure 3 .

 3 Figure 3.7A for each chromosome separately (Appendix Figure B.1). The results obtained for the entire genome were recapitulated for each chromosome, with minor exceptions detailed in the legend of Figure B.1. To quantify these results, we computed the correlation between the genome wide C RF D 10kb and the correlation matrix of RFD for each chromosome C chr i RF D 10kb using the procedure dened in section 2.3.3 . The correlation between the global RFD correlation matrix and each chromosome's RFD correlation matrix was high (>0.893), albeit sensitive to chromosome size as expected (rst raw of matrix and last column in the matrix of Figure 3.8). We similarly computed the pairwise correlation coecient between individual chromosomes RFD correlation matrices (Figure 3.8). Correlations between individual chromosomes were also high (>0.79), but again sensitive to chromosome size (Figure 3.8).

Figure 3 . 8 :

 38 Figure 3.8: Replication changes between cell lines are widespread through the genome. Correlation coecients between the global genome correlation matrix between RFD proles shown in Figure 3.7 A and the individual chromosomes' correlation matrices shown in Figure B.1. Pearson correlation coecient values are color-coded from blue (0.4) to red (1) using the color bar on the right.

Figure 3 . 9 :

 39 Figure 3.9: RFD proles are more conserved in high GC-content regions.(A-E) Correlation matrix of RFD proles depending on the GC content; 10 kb windows were grouped in GC-content categories following the 5 isochores classication of the human genome in light isochores L1 (GC ≤ 37 ; C L 1 RF D ; A) and L2 (37 ≤ GC < 41; C L 2 RF D ; B), and heavy isochores H1 ( 41 ≤ GC < 46; C H 1 RF D ; C), H2 (46 ≤ GC < 53 ; C H 2 RF D ; D) and H3 (GC ≥ 37 ; C H 3 RF D ; E); coverage (Cov) of the sequenced human genome by each isochore family is provided in the header of each column; Pearson correlation coecient values are colour-coded from blue (0.4) to red (1) using the colour bar on the right (Materials and Methods).

Figure 3 . 10 :

 310 Figure 3.10: Random 5 Mb probes are sucient to recover the cell line classication. (Blue, boxplot) Distribution of the correlation coecient between the RFD correlation matrix for the complete genome (Figure 3.7A) and 1000 random probes of indicated of scale (Mb) consisting of 5 to 5,000 randomly located 10 kb windows. (Orange), probability of observing the same 3 groups classication of cell lines as in Figure 3.7A, when performing the classication using a random mode of the indicated scales.

Figure 3 .

 3 Figure 3.11: RFD proles are more conserved in high GC-content regions. (A-E) Correlation matrix of RFD proles depending on the GC content; 10 kb windows were grouped in GC-content categories following the 5 isochores classication of the human genome in light isochores L1 (GC ≤ 37 ; C L 1 RF D ; A) and L2 (37 ≤ GC < 41; C L 2 RF D ; B), and heavy isochores H1 (41 ≤ GC < 46; C H 1 RF D ; C), H2 (46 ≤ GC < 53 ; C H 2 RF D ; D) and H3 (GC ≥ 37 ; C H 3 RF D ; E). (F-J) Matrices of correlation dierences ∆C I RF D =C I RF D -CRF D where I can be L1, L2, H1, H2 or H3; correlation dierence values are color coded blue (resp. red) for negative (resp. positive) dierences using the color bar on the right; a blue (resp. red) color indicates that RFD proles are less (resp. more) conserved in the considered isochore than in the 22 autosomes. Matrices row and column order is the same as in Figure 3.7A.

Figure 3 . 12 :

 312 Figure 3.12: Analyses of the scale dependence of the mean correlation dierence between isochore-and global genome-based RFD correlation matrices. The mean values of correlation dierences ∆C I RF D = C I RF D -CRF D(where I can be L1, L2, H1, H2 or H3) were computed for each isochore as in Figure3.11. F-J at 10kb, 100 kb, 200 kb and 1 Mb scales. The mean correlation dierences increase with GC-content at all analysis scales.

The 5

 5 Mb scale was chosen because at this scale statistical noise is small enough to not interfere with cell lines classication, based on the result of the random probes approach (Figure 3.10). We computed for each non overlapping 5 Mb window the RFD correlation matrix. We observed, for example, among all 5 Mb non overlapping windows of chromosome 1 (Appendix, Figures A.4) that the classication of the 12 cell lines using RFD proles is not conserved in the majority of windows (totally blue or red matrices). Even though some 5 Mb regions that conserved the classication of cell lines in three groups. Moreover, we observed a strong association within the myeloid cell lines along the genome, illustrated by a red square in the majority of the RFD correlation matrices (Appendix, Figure A.4). In Figure 3.13A, we represented a 5 Mb window selected from the chromosome 1 (35-40 Mb) where the RFD proles of 12 cell lines are similar. The correlation matrix of these RFD proles (Figure 3.14B) illustrates that the RFD proles among 12 cell lines are very correlated and that cell line classication according to 3 types is lost due to replication program proximity. In contrast, a disturbed RFD correlation matrix corresponding to variable RFD proles among the 12 cell lines is observed in the region between 185 Mb and 190 Mb of chromosome 1 (Figure 3.14).In Figure3.15 we compare the correlation coecient between global correlation matrix and RFD correlation matrix when using windows rather than probes (3.10) (5 kb to 50 Mb).For all scales between 5 kb and 50 Mb, we observe that the correlation of the complete genome RFD correlation matrix with the RFD correlation matrix computed in windows is smaller than with RFD correlations obtained with random probes of the same sizes. For example, for the distributions of the correlation coecients of 5 Mb windows (Green boxplot) and 5 Mb probes (Blue boxplot) we see that the distribution of correlation coecients of the rst one vary between ∼0.5 and ∼0.75 and vary between ∼0.95 and ∼1 in the second one (Blue boxplot) (Figure3.15). This suggests that we have more heterogeneity of RFD proles changes among the 12 cell lines in 5 Mb window than expected scale 5 Mb but also at for all scales >50 Mb.

Figure 3 . 13 :AFigure 3 . 14 :Figure 3 . 15 :

 313314315 Figure 3.13: Example of Mb window with a highly stable replication program. (A) RFD of the 12 cell lines in a 5 Mb regions. (B) Correlation matrix of the RFD proles presented in A.

Figure 3 . 16 :

 316 Figure 3.16: Distribution of Z-score values. We represented the histogram of the ratio Mean by standard deviation of the dierences correlation matrix between each 5Mb window and the global one. Blue hatched the 5% unstable regions, red hatched the 5% stable regions of genome.

Figure 3 .Figure 3 . 17 :

 3317 Figure 3.16 showed that we have a large distribution of Z-scores, thus in many 5 Mb window regions the global cell lines classication based on RFD proles is disturbed. To validate this method, we computed the Z-score for the GC-rich (H3, 4% of genome, Figure3.11J) regions where the RFD proles among the 12 cell lines were estimated as stable and we obtained 1.76, which corresponds to a stable RFD proles among the 12 cell lines. In contrast we computed the Z-score for the GC-poor regions (L1, 26% of genome, Figure3.11F) and we obtained -1.86, which corresponds to variable RFD regions in agreement with the previous results (section 3.4.1). Then, we selected the very similar regions by taking the highest 5% Z-score values (z>1.64). We selected the 5% lowest Z-score (z< -1.25) to identify the regions where the RFD proles vary widely among the 12 cell lines (Figure3.16). Then, to obtain further characterization of these regions, we compared the GC, FPKM 6 and MRT between the detected regions of stable and variable replication program (Figure3.17). The regions with6 Fragments Per Kilobase of exon model per Million mapped fragments is a unit of gene expression level as estimated by RNA-seq proling

Figure 3 . 18 :

 318 Figure 3.18: MCR distribution. (A) Cumulative distribution function of MCR computed in non overlapping 500 kb windows (Eq 3.4). (B) MCR proles along the exemplary 5 Mb regions with a stable or variable replication program. (Top) Same stable regions as in Figure 3.13. (Bottom) Same variable regions as in Figure 3.14. MCR prole was computed in 500 kb windows centered every 10 kb. Color coding dened in A.

  .18A. We partitioned the genome in 4 classes of equal size based on the MCR values: very variable regions for MCR<0.46, moderately variable regions for MCR ∈ [0.46, 0.52[, moderately stable regions for MCR ∈ [0.52, 0.57[ and very stable regions for MCR≥0.57. To validate the MCR score, we computed the global mean dierence of RFD proles matrix among the 12 cell lines at 10 kb according to the MCR level. The dierence between each mean dierence RFD proles matrix according to each MCR level and the global mean dierence RFD proles matrix among the 12 cell lines show an evolution of color from blue (low MCR) with dierences larger than average to red (high MCR) with dierences smaller than average (Figure 3.19). Thus, we conrmed that the regions with high MCR are the regions with less RFD dierences. In contrast, regions with low MCR are regions with high RFD dierences among the 12 cell lines.

Figure 3 . 19 :

 319 Figure 3.19: Evolution of RFD correlation matrix with MCR score. 10 kb windows were grouped in 10 categories of equal size based on the MCR value of the 500 kb regions centered on the windows. (A-J) Correlation matrix of RFD proles depending on the MCR level. (A'-J') Matrices of dierence of correlation matrix. Matrices row and column order is the same as in Figure 3.7A.

Figure 3 . 21 :

 321 Figure 3.21: High correlation between MCR level and GC content.. Cdf of the GC content at 100 Kb according to the MCR level at 500 kb. From green, variable RFD proles regions, to red highly stable RFD proles regions.

Figure 3 . 22 :Figure 3 . 23 :

 322323 Figure 3.22: The global MCR score are highly correlated to pairwises cell type MCR scores.. (ABC) 2D probability distribution function (pdf) between the global MCR score and MCR score computed in pair of cell lines type M CRML (myeloid and lymphoid), M CRMA (myeloid and adherent), M CRLA (lymphoid and adherent). (DEF) 2D pdf between the MCR computed per cell lines type ( M CRM for myeloid cell lines, M CRL for lymphoid cell lines, M CRA for adhrent cell lines) and the global MCR score . Color-bar corresponding the the probability density function.

Figure 3 . 24 :

 324 Figure 3.24: Distribution of local extremum values of RFD slopes proles in each cell lines. Pdf of the slopes at the local maximum of RFD slopes proles (blue). Pdf of the opposites of the slopes at the local minimum of RFD slopes proles (green).

3. 5

 5 Preferential IZ and DNA replication program 3.5.1 Recovering preferential replication initiation zones IZ are ascending segments of the RFD proles (Figures 1.7 and 3.4, [3]) and can thus be detected as regions of maximal values of the RFD proles slopes. A concern about experimental

Figure 3 . 25 :

 325 Figure 3.25: Number of common IZ with slope>1%RFD per kb at 30 kb resolution for all pairs of cell lines.The upper triangular part of the matrix represents the observed count of common origin with slope>1%RFD per kb between pair of cell lines in 30 kb windows. The lower triangular part of the matrix represents the expected count of common IZ computed >1%RFD per kb in each cell lines. Main diagonal is the number of detected IZ of slope>1%RFD per kb in each cell lines.

Figure 3 . 26 :

 326 Figure 3.26: Matrix of observed/expected numbers of common for each cell line pair for a threshold of 1%RFD per kb. Each square of this matrix represents the normalized count of common origin between pair of cell lines. Normalization was done by dividing the observed count of common IZ (upper part of the matrix in Figure A.7) by the expected one (lower part of the matrix in Figure A.7).

  between 1943 common IZ (between TLSE19 and IARC385) and 4226 (between TF1_GFP and TF1_BcrAbl_1M) (Figure A.7). The low number of common IZ between TLSE19 and IARC385 corresponds the low RFD proles correlation coecient between them as shown in the Figure 3.7A. The high number of common IZ detected between TF1_GFP and TF1_BcrAbl_1M is associated with the high RFD correlation coecient between them. In addition, we see that the number of common IZ between the TF1 cell line and the K562 also increases as expected from CML progression (CML progression system will be specically analyzed in Chapter 5). To validate the detected numbers we computed the expected number of common IZ between each pair of cell line at 30 kb windows as the mean for random distribution of IZ using Binomial law (Figure A.7). About 300 common IZ are expected between each two cell lines (Figure A.7) so is smaller than the observed counts. In Figure 3.26, we represented the ratio between the observed and expected common IZ between all pairs of cell lines computed by dividing the upper triangular part of the matrix (Figure A.7, observed) by the lower triangular part of the matrix (Figure A.7, expected). We systematically observed in average about >6 times more common IZ than expected. Except for HeLa and for IB118 cell line, the global matrix structure showed in the Figure 3.26 is very similar to the structure of the RFD correlation matrix showed in the previous section of this chapter (Figure 3.8 A). Thus, in each type of cell lines, the cell lines share with high concordance the same IZ. Moreover, we can see the relation between the K562 and the CML progression as in the global RFD correlation matrix (Figure 3.26). For the lymphoid cell lines, we observe a highly conserved IZs between the two BL cell lines (BL79 and Raji), furthermore IARC385 shares more IZ with BL cell lines (BL79 and Raji) than GM06990. Similar results were obtained when increasing the threshold to 2 (Figure A.[START_REF] Alberts | Essential cell biology[END_REF]. We observed between ∼16

3 .Figure 3 . 27 :

 3327 Figure 3.27: High density of replication origins in high GC content replication stable regions. (A) The IZ density detected using the local maximum of RFD slope (MS) proles in 5 categories of GC at 100 kb from poor to rich ([0,0.37], [0.37,0.42], [0.42,0.47], [0.47,0.53],[0.53,1]). Blue dots all MS>0.5%RFD per kb (S1), orange dots MS>1%RFD per kb (S2) and green dots the MS>2%RFD per kb (S3). (B) Cdf of GC at 100 kb in each group of slope s(Blue line for S1, Orange line for S2, Green line for S3, Black line in the whole genome). (C) Same gure as in A but for MCR at 500 kb from variable to conserved RFD proles. (D) Cdf of MCR at 500 kb in each group of RFD slope local maximum.

1 -

 1 light R bands[START_REF] Holmquist | Characterization of giemsa dark-and light-band DNA[END_REF] associated to GC-rich and highly transcribed regions, the second half corresponding to the Giemsa-dark G bands associated to GC-poor and low transcriptional activity regions[163165]. In the 2000s, the emergence of microarray technologies allowing genome wide DNA replication timing analyses in many organisms as drosophila[START_REF] Schübeler | Genome-wide DNA replication prole for drosophila melanogaster: a link between transcription and replication timing[END_REF][START_REF] Macalpine | Coordination of replication and transcription along a drosophila chromosome[END_REF] and human[START_REF] White | DNA replication-timing analysis of human chromosome 22 at high resolution and dierent developmental states[END_REF][START_REF] Woodne | Replication timing of the human genome[END_REF][START_REF] Jeon | Temporal prole of replication of human chromosomes[END_REF]] conrmed these correlation. We conrmed these observations at the chromosomal scale. We observed a strong correlation of the chromosomes average MRT with their GC content, their average FPKM and DNase 2 sensitivity (Figure4.1A). Note that, the correlation coecient between MRT and DNase sensitivity (r=-0.605) was stronger than the correlation between MRT and FPKM (r=-0.56) (Figure4.1B).Generally, active and open euchromatin regions are replicated early in S phase, and repressed and closed heterochromatin regions are replicated late in S phase[START_REF] Rhind | DNA replication timing[END_REF]. In the mammalian genome the replication timing is also associated with the subnuclear position (Figure1.5). Using logistic regression[START_REF] Woodne | Replication timing of the human genome[END_REF] [START_REF] Woodne | Replication timing of the human genome[END_REF] found a signicant positive correlation

Figure 4 . 1 :

 41 Figure 4.1: Association of Mean Replication Timing with GC content, transcription and Dnase sensitivity in K562. Mean of GC content (A) FPKM (B), and DNase sensitivity (C) of each chromosome as a function of the mean of MRT in each chromosome. The orange line represent the linear regression. The r value in the top of each panel represents the correlation coecient of Pearson.

Figure 4 . 2 :

 42 Figure 4.2: Genome-wide proling of replication timing. Mean Replication Timing (MRT) along a 20 Mb segment of chromosome 3 in 2 lymphoid (top), 3 myeloid (middle) and 2 adherent (bottom) cell lines, as indicated in the top right corners of each panel.

Figure 4 . 3 :

 43 Figure 4.3: Genome-wide proling of gene expression. FPKM at 200 kb (Eq. 2.2) along a 20 Mb segment of chromosome 3 in 4 lymphoid (top), 4 myeloid (middle) and 4 adherent (bottom) cell lines, as indicated in the top right corners of each panel.

Figure 4 . 4 :

 44 Figure 4.4: RNA-seq biological replicate classication based on correlations between gene expression computed at 200 kb (Figure 4.3). (Bottom) Correlation matrices between RNA-seq proles. Pearson correlation coecient values are color-coded from blue (0.4) to red (1) using the color bar on the right. (Top) Dendrogram representation of the hierarchical classication of RNA-seq experiments based on the corresponding correlation matrix; ordinate is the correlation distance.

Figure 4 . 5 :

 45 Figure 4.5: Cell line classication based on correlations between gene expression proles. (Bottom) Correlation matrices between RNA-seq. CF P KM Pearson correlation coecient values are colour-coded from blue (0.4) to red (1) using the color bar on the right (Materials and Methods). (Top) Corresponding dendrogram representation of the hierarchical classication of cell lines; ordinate is the correlation distance (Eq. 2.7)

  , 4.5). By RNA-seq, the two LMS IB118 and TLSE19 were more correlated to each other than to IMR90 and less correlated to HeLa, which in fact clustered with myeloid cells. The cells of origin and driver mutations of LMSs are currently unclear. These results

Figure 4 . 6 :

 46 Figure 4.6: Cell line classication based on the MRT proles. (A) Correlation matrices between MRT proles ((CMRT ; A); Pearson correlation coecient values are color-coded from blue (0.4) to red (1) using the color bar on the right (Materials and Methods). (Top) A corresponding dendrogram representation of the hierarchical classication of cell lines is shown on top of each correlation matrix; ordinate is the correlation distance. (B) Cumulative distributions of the absolute MRT changes (|∆MRT|) between cell lines. Each curve is color-coded according to the pair of cell lines indicated in the insert. The considered threshold of signicance [111] (|∆MRT| > 0.2) is indicated by a vertical dotted line.

  and 4.6 for each chromosome separately (Figure B.1).

Figure 4 . 7 :

 47 Figure 4.7: Transcription program conservation for regions with dierent GC-content. (A-E) Correlation matrix of RNA-seq proles depending on the GC content; 10 kb windows were grouped in GC-content categories following the 5 isochores classication of the human genome in light isochores L1 (GC ≤ 37 ; C L 1 RN A-seq ; A) and L2 (37 ≤ GC < 41; C L 2 RN A-seq ; B), and heavy isochores H1 (41 ≤ GC < 46; C H 1 RN A-seq ; C), H2 (46 ≤ GC < 53 ;C H 2 RN A-seq ; D) and H3 (GC ≥ 37 ; C H 3 RN A-seq ; E). (F-J) Matrices of correlation dierences where can be L1, L2, H1, H2 or H3; correlation dierence values are color coded blue (resp. red) for negative (resp. positive) dierences using the color bar on the right; a blue (resp. red) color indicates that RFD proles are less (resp. more) conserved in the considered isochore than in the 22 autosomes. Proportion of genes (Gene prop.) within each isochore family is provided in the header of each column

Figure 4 . 8 :

 48 Figure 4.8: MCR in late and early timing regions. The lines represent the cdf of MCR in consistently late replicating regions ( M RT >0.7), the dashed lines represent the cdf of MCR in consistently early replicating regions ( M RT <0.3). Each color correspond to a scale of computation of MCR , as indicated.

Figure 4 . 9 :Figure 4 . 10 :

 49410 Figure 4.9: Identication of large (>1.5 Mb) stable DNA replication timing program domains. We represent the RFD proles at 10 kb (C, resp. D), MCR at 500 kb and Mean Replication Timing (MRT) at 100 kb (A, resp. B) of 7 (resp. 12) Mb of early (resp. late) replicating regions of chromosome 2 (resp. chromosome 4). MCR code is dened in Figure 3.18. Other examples of large early and late replication domains are provided in Figure B.4.

Figure 4 .

 4 Figure 4.9A (resp. B) large domains (>1.5 Mb) of consistently early (resp. late) MRT. The large early replicating domain (>1.5 Mb) is associated to a high MCR values (Figure 4.9A) and the large late replicating domains is associated to a low MCR values (Figure 4.9B). Thus we observed an association between the large replication timing domains and large domains (>1.5 Mb) of constant MCR. To investigate if this nding is a global result, we associated M RT at 100 kb to all high MCR (MCR>0.55) and low MCR (MCR<0.47) proles at 500 kb (Figure 4.10C, blue). Then we associated the M RT to large domains of high and low MCR (>1.5 Mb) (Figure 4.10C, orange). We observed that the large domains of low MCR are concentrated in late replicating regions ( M RT >0.6, orange dot). But the M RT in large domains of high MCR are more dispersive, so we can have some regions which are stable in terms of RFD proles but are replicating late. Thus, the large domains of variable RFD proles are majoratively associated to late replicating regions. In Figure 4.10A, we selected all large domains (>1.5 Mb) of late and early replication and we associated them to MCR values at 500 kb. We observed that the late replicating domains are associated to low MCR values (orange distribution), in opposite the early replicating domains are associated to high MCR values (blue distribution). The reciprocal analysis is shown in the Figure 4.10B, we associated the M RT to the large domains (>1.5 Mb) of variable (MCR<0.47) and stable (MCR>0.55)RFD proles. We conrmed the very strong association between the variable RFD proles domains and the high M RT values (orange distribution). The stable RFD proles domains are mostly associated to early replicating regions. However, the latter association is weaker than the former. These results suggested that globally we have a strong association between large domains (> 1.5 Mb) of late (resp. early) replication timing and large domains (> 1.5 Mb) of variable (resp. variable) RFD proles. Finally, in Figure4.9CD, we represented the RFD proles at 10 kb for three cell lines (GM06990, TF1_GFP and IMR90) that represent each type of cell line group, in respectively early and late constant timing regions. The large

Figure 4 .Figure 4 . 11 :Figure 4 . 12 :

 4411412 Figure 4.11: MRT, MCR and RFD slope proles in 10 Mb early stable replicating regions. (A) MRT proles of 5 cell lines at 100 kb in 10 Mb early replicating regions. (B) MCR proles at 500 kb in stable RFD proles. (C) RFD slope proles for the 12 cell lines.

Figure 4 . 13 :

 413 Figure 4.13: Only strong IZ in GM06990 are conserved in Raji. (A) (resp. B) Scatter plot of the RFD slopes in Raji associated to the local maximum of RFD Slope (SM) in GM06990 >1.Blue dots represent the early replicating regions according to the MRT of GM06990 at 100 kb (resp. stable replicating regions according to MCR at 500 kb). Orange dots represent the late replicating regions according to the MRT of GM06990 at 100 kb (resp. variable replicating regions according to MCR at 500 kb). Grey dots represent the rest of the genome.

Figure 4 . 14 :

 414 Figure 4.14: High density of origin in early replication program regions. (left) Replicative density of IZ detected using the slope of RFD proles in 5 categories of MRTTF1_GFP at 100 kb from early to late ([0, 0.2], [0.2, 0.4], [0.4, 0.6], [0.6, 0.8], [0.8, 1]). Blue dots for all slopes>0.5%RFD per kb (S1), orange dots slopes>1%RFD per kb (S2) and green dots for slopes>2%RFD per kb (S3). (right) Cdf of MRTTF1_GFP in each group of slope (Blue line for S1, Orange line for S2, Green line for S3, Black line for MRTTF1_GFP at 100 kb in the whole genome).

Figure 4 . 15 :

 415 Figure 4.15: Transcription associated to replication change. (A) We represent the cumulative distribution function of the average FPKM value on the 12 cell lines according to their MCR level. (B) We represent the cdf of the relative changes between the pair of cell line with the lowest correlation in the RFD proles (which is IARC385 and TLSE19, Figure 3.7A).

Figure B. 5 )Figure 4 . 16 :

 5416 Figure B.5). This Figure B.5 conrms that we have a global association between the replication program changes and the relative transcription program changes but in opposite directions.

Figure 4 . 17 :

 417 Figure 4.17: Average absolute dierence matrices of FPKM proles depending on the MCR level. 10 kb windows were grouped in MCR categories following the 4 groups classication of the human genome in extremely stable RFD proles regions with MCR≥0.57 (G1, A) , the moderately stable RFD proles regions with 0.52≤MCR<0.57 (G2, B), the moderately variable RFD proles regions with 0.42≤MCR<0.52 (G3, C) and the extremely variable RFD proles regions with MCR<0.42 (G4,D); (F-H) Matrices of dierences between average absolute FPKM dierence matrices in each MCR levels and the genome-wide matrix (Figure 4.16)

Figure 4 . 18 :

 418 Figure 4.18: Average relative dierences between gene expression proles. Same as 4.16 but for relative dierences.

Figure 4 . 19 :

 419 Figure 4.19: Average relative dierence matrices of FPKM proles depending on the MCR level. Same as 4.17 but for relative dierences.

  in GM06990 and Raji (1753) Not Expressed in GM06990 (25) Not Expressed in Raji (265) All promoter

Figure

  Figure 4.20: Association between gene expression and RFD proles in GM06990 and Raji. (A) RFD proles and (B) FPKM proles of GM06990 (blue) and Raji (green) computed in 10 kb windows along a 5 Mb region of chromosome 2. (C) Cdf of the MCR values at the dierentially expressed gene promoters (red) and at all promoters (black). (D) Cdf of the length of gene expressed in both GM06990 and Raji (orange), all genes (black), dierentially expressed genes (red), expressed in Raji only (green) and expressed in GM06990 only (silver). (E) Counts of MCR level at the promoters of gene expressed in GM06990 and Raji (orange), all promoters (black), dierentially expressed genes (red), genes expressed in Raji only (green) and genes expressed in GM06990 only (silver).

  Figure B.7) showed that the HeLa cell line is close to the adherent cell line group and myeloid cell line group. This is very compatible with the RNA-seq correlation matrix (Figure 4.5) where HeLa was classied with myeloid cell line instead of Adherent cell lines.

Figure 4 .

 4 20E). More importantly, we observed that dierentially expressed genes have a lower MCR at their promoters compared to the complete gene set (Figure4.20C). Therefore, the signicant gene expression changes are associated with late replicating regions based on the association between the large domains of late replication timing values and low MCR domains (section 4.3.1). This conrmed our observation that low MCR regions are associated to the largest relative changes of expression (section 4.4.1).

Figure 4 . 21 :

 421 Figure 4.21: Association between changes in replication initiation eciency and in gene expression. (A) The slope dierences between GM06990 and Raji (∆slope Raji GM 06990 = slopeGM06990-slopeRaji) at the promoter of signicantly expressed genes with FPKMGM06990 > FPKMRaji (green), with FPKMGM06990 < FPKMRaji (black), with no signicant gene expression change in green. (B) Same analysis between TF1_GFP and TF1_BcrAbl_1M.

Figure 4 . 22 :

 422 Figure 4.22: Type of replication changes. We illustrated the 3 types of RFD proles dierences between GM06990 (blue) and Raji (Orange). RFD proles were computed at 10 kb. The red arrows correspond to the orientation of the genes. From left to right : (Type I) corresponds to GM06990 specic IZ, (Type II) corresponds to Raji specic IZ and (Type III) corresponds to the common IZ detected in the two cell lines that may change eciency or shift position

Figure 5 . 1 :

 51 Figure 5.1: Early CML progression model. The CML, progression model includes the TF1_GFP cell line, the TF1_BCRABL_1M cultured 1 month after the transduction of the BCR-Abl1 gene and TF1_BCRABL_6M cultured 6 months after the transduction of the BCR-Abl1 gene. Step 1 (resp. 2) changes correspond to the changes of RFD proles between TF1_GFP and TF1_BcrAbl_1M (resp. TF1_BcrAbl_1M and TF1_BcrAbl_6M)

Figure 5 . 2 :

 52 Figure 5.2:The largest RFD changes induced by 1 month of BCR_ABL1 expression are observed in GC-poor, lowly-expressed and late replicating regions. Cumulative distribution functions (cdf) of GC content (A), transcription in TF1_BrcAbl_1M F P KM1M (B) and MRT in TF1_BrcAbl_1M (C) computed in non-overlapping 200 kb windows of the 22 autosomes. Cdfs were determined for all windows (all, black) or selected windows with the 20% (blue), 10% (green), 5% (yellow) and 1% (red) largest RFD changes between TF1_GFP and TF1_BrcAbl_1M.

Figure 5 . 3 :

 53 Figure 5.3: The largest RFD changes between 1 month and 6 months of BCR_ABL1 expression are observed in GC-poor, lowly-expressed and late replicating regions. Cumulative distribution functions (cdf) of GC content (A), transcription in TF1_BrcAbl_1M (B) and MRT in TF1_BrcAbl_1M (C) computed in non-overlapping 200 kb windows of the 22 autosomes. Cdfs were determined for all windows (all, black) or selected windows with the 20% (blue), 10% (green), 5% (yellow) and 1% (red) largest RFD changes between TF1_BrcAbl_1M and TF1_BrcAbl_6M.

Figure 5 . 4 :

 54 Figure 5.4: The cell lines clustered in accordance to CML progression.(A-C) Correlation matrices between RFD proles in accordance to CML progression ( A), RNA-seq (B) and MRT proles (C); Pearson correlation coecient values are colour-coded from blue (0.4) to red (1) using the colour bar on the right. (Top) A corresponding dendrogram representation of the hierarchical classication of cell lines is shown on top of each correlation matrix; ordinate is the correlation distance.

Figure 5 . 5 :

 55 Figure 5.5: Manual annotation of RFD prole changes during the rst two steps of the CML progression model.(A-D) RFD proles computed in 10 kb non-overlapping windows for TF1_GFP (black), TF1_BrcAbl_1M (yellow) and TF1_BrcAbl_6M (blue) are visualized in 2 Mb regions along the 22 human autosomes. (A) Region along chromosome 12 where RFD proles in the 3 cell lines do not present any signicant dierence. (B) Region along chromosome 2 illustrating IZs whose eciency is repeatedly enhanced (green) or weakened (purple) in Steps 1 and 2 of the model CML progression. (C) Region along chromosome 18 having two IZ eciency changes in Step 1 (enhanced, leftmost green; weakened, purple) that are conrmed in Step 2. The rightmost green line marks an IZ which is repeatedly enhanced at Steps 1 and 2 (as in B). (D) Region along chromosome 5 illustrating an IZ activated during Step 1 (red line) and silenced during Step 2.

5. 4

 4 Manual annotation of initiation zones eciency changes RFD proles changes in the rst two steps of the CML initiation and progression model (Step 1: TF1_GFP -> TF1_BCRABL_1M (Figure 5.1); Step 2: TF1_BCRABL_1M -> TF1_BCR_ABL_6M (Figure 5.1)) were annotated by manually scanning the proles in 2 Mb windows. An IZ present in the initial state was annotated as Silenced if no IZ was present in the following state. Weakened when IZ eciency decreased between the two consecutive states (Figure 5.5C). Enhanced when IZ eciency increased between the two consecutive states (Figure 5.5B). IZs inactive in the initial state but active in the nal state were annotated New (Figure 5.5D). In contrast, we do not observed any change in MRT proles (Appendix V Figure C.2A) or transcription prole (Appendix V Figure C.2B). The regions where we have a change are regions desert of active genes and late replicating. Moreover, each locus annotated for a Step 1 or Step 2 change (total 1027) was also annotated for its status in Step 3 (TF1_BCRABL_6M -> K562). Step 3-specic changes were too numerous to be manually annotated. The 1027 manually annotated loci included 253 IZs which changed eciency at Step 1 but not Step 2, 551 IZs eciency changes during Step 2 but not Step 1, and 223 IZs eciency changes in both Step 1 and Step 2. In total, this database encompasses 476 and 774 and 716 eciency changes in Steps 1 and 2 and 3, respectively (Table5.1).5.5 BCRABL1 expression continuously induces replication changes in gene desert5.5.1 The targeted initiation zones were more frequently weakened than enhanced We focused further our analysis on the early CML progression model, where changes in RFD can be directly attributed to changes in expression of the BCR_ABL1 oncogene in TF1 cells. Consistent with their high global correlation coecients (>0.95), we observed a striking identity of the RFD proles of the three cell lines over most of the genome, as exemplied in Figure 5.5A. This facilitated the detection and manual annotation of IZ eciency changes, scored as new, enhanced, weakened or silenced IZs, at each step of CML progression (section 5.4). A few examples of annotated 2Mb segments are shown on Figure 5.5. RNA-seq and MRT proles of the same regions are shown on Figure C.2. Along the genome, 476 changes during Step 1 and 774 changes during Step 2 were observed (Table 5.1). The distributions of IZ eciency changes were strikingly similar at Step 1 and Step 2 (Figure 5.6A). Weakened IZs were by far the most frequent at each step (55% and 66%, respectively) and enhanced IZs the second most frequent. We then analyzed the evolution of the 476 Step 1 changes during Step 2 (Figure 5.6B). The behavior of the Step 1 changes during Step 2 signicantly depended on the type of change at Step 1 (P < 0.001 using a χ 2 test of independence

Figure 5 . 7 :

 57 Figure 5.7: Weakened initiation zones eciency changes in step 2 are more associated to conserved RFD proles. Cumulative distribution function of MCR at 500 Kb in accordance to the IZ changes in step 1 (left) and step 2 (right). Red lines for New ring origin. Blue lines for silenced origin. Purple for Enhanced IZ, and Green lines for Weakened IZ.

Figure 5 . 8 :Figure 5 . 9 : 2 (

 58592 Figure 5.8: Initiation zones eciency changes in response to 1 month of BCR_ABL1 expression are observed in GC-poor, lowly-transcribed and late replicating regions. Cumulative distribution functions (cdf) of GC content (A), transcription in TF1_BCRABL_1M (B) and MRT in TF1_BCRABL_1M (C) computed in non-overlapping 200 kb windows of the 22 autosomes. Cdfs were determined for all windows (all, black) or limited to windows with Silenced (blue), Weakened (green), Enhanced (violet) and New (red) IZ. Similar results are obtained using transcription and MRT data of TF1_GFP (Figure C.3).

  .8 we associated for the step 1 changes categories (New, silenced, Enhanced and weakened) the GC content, MRT and genes expression level represented by the FPKM of TF1_BcrAbl_1M at 200 kb. We conrmed that Step 1 changes preferentially occurred in GCpoor, lowly expressed and late replicating regions (Figure5.8). Interestingly, this tendency was more pronounced for new or silenced IZs than for weakened or enhanced IZs, in other words for more extreme changes. Similar results were obtained whether RNA-seq and MRT data from TF1_BCRABL_1M (Figure5.8) or from TF1_GFP (Appendix, Figure C.3) were used. A more complex situation was observed for Step 2 changes (Figure 5.9). New, enhanced and silenced IZs again concentrated in GC-poor, lowly expressed and late replicating regions, but weakened IZs now were more often found in GC-rich, highly expressed, early replicating DNA, coherently with their preferential location in high MCR regions at this step (Figure 5.7). To address whether weakened IZs at Step 2 were associated with nearby gene transcription changes, we plotted the RNA expression (by 200 kb windows) ratio in TF1_BCRABL_6M over TF1_BCRABL_1M, as a function of their mean expression level (Figure 5.10A, MA plot) and computed the cumulative distribution of the expression changes (Figure 5.10B), for the total genome or for windows containing at least one weakened IZ at Step 2. The results indicate that IZ weakening events at Step 2 were signicantly associated with nearby transcription repression. We observed that correlation coecients between cell lines were generally smaller by RFD (Figure 3.7A) than by RNA-seq (Figure 4.5). We investigated in the CML system whether, when focusing on early replicating regions, changes in RNA predicted changes in RFD and IZs. The distribution of RNA expression changes in early replicating windows with the largest RFD changes at Step 2, was shifted towards RNA repression when compared to all early windows (Figure 5.11AB), consistent with the association of weakened IZs with transcription repression (Figure 5.10). However, we found no dependence of the distribution of Step 2 RFD changes on RNA expression changes in early replicating regions (Figure 5.11C).

log 10 (Figure 5 . 10 :

 10510 Figure 5.10: Weakened initiation zones between 1 month and 6 months of BCR_ABL1 expression are associated with transcription repression.(A) log10-ratio of FPKM values in TF1_BCRABL_6M over TF1_BCRABL_1M as a function of their log10 (geometric) mean. FPKM values were computed in non-overlapping 200 kb windows. Only windows expressed in both cell lines were considered (FPKM > 0.01). Windows containing at least one weakened IZ at Step 2 are in orange, the rest is in blue. (B) Cdf of the FPKM log10-ratios for windows with FPKM (geometric) mean > 0. Blue and orange, as in (A).

3. 7 )

 7 and MRT (Figure 4.6) proles are clustered in three separate groups corresponding to lymphoid, myeloid and adherent cells. Similar results were obtained by RNA-seq (Figure 4.5), except that HeLa clustered with myeloid instead of adherent cells.

Figure 4 . 13 .

 413 Figure 4.13. This suggested a cell line specic control of IZ in late replicating domains during normal and pathological dierentiation. Late replication are associated to low gene density suggesting a decoupling between replication program changes and gene expression changes, at least in these regions.

(Figure 3 . 3

 33 around position 5 Mb). A possible explanation of this observation is that OF biochemistry could depend on genomic or epigenomic context. This prompt our experimental collaborators at IBENS, Paris to extend the OK-seq protocol to paired-end sequencing. Indeed, these new Ok-seq experiments should in principle enable us to reconstruct the OF size statistics at every locus as the distribution of paired-end reads outer length (Figure D.1). Paired-end Ok-seq experiments performed on previously proled cell lines (GM06990, Raji, BL79 and IARC385) already demonstrated that this now approach is highly compatible with the original protocol (Figure D.5) and that at least for 3 cell lines (Raji, BL79, IARC385) the observed genome-wide estimated OF length distribution is compatible with the expected OF full size of ∼ 150 bp (Figures D.2, D.3 and D.4). This new data demonstrate the potential of Ok-seq for genomic analyses of the dynamics of OF processing. OK-seq experiment can also be designed to ask specic question about the implication of oncogenes on RS. The Myc expression in hematopoietic cells can induce tumorigenesis. Myc over-expression is believed to be one of the primary events in the malignant transformation of Burkitt's lymphoma. Overexpression of Myc increases cell size and energy production, induces proliferation and genomic instability, and favors apoptosis. Moreover, Myc over expression can cause RS from the earliest tumorigenesis. New data constructed to question the role of C-myc overexpression on RS are new available (Hyrien team, IBENS, Paris). They are based on a B lymphocyte cell line (P493) where EBNA and Myc where transfected under the control of inducible promoters (Figure D.6). EBNA is required to assure cell proliferation in the absence of Myc expression. The P493 cell line was obtained from a B-cell line by allowing expression of one or both of the EBNA and Myc (Table D.1) and each culture was proled by Ok-seq. The RFD proles of P493 cultures with the expression of EBNA and Myc or with EBNA expression alone are very correlated to each other but are less correlated to the RFD proles with Myc expression alone. Hence, EBNA expression appears to be a stronger factor inuencing RFD proles than Myc (Figure D.[START_REF] Watson | Genetical implications of the structure of deoxyribonucleic acid[END_REF]. This shows that OK-seq proling has the potential to inform us on the eect of Myc on the replication program, even though in the current system there will be a challenge to factor out the eect of EBNA expression.To conclude, in Chapter 4 we have identied 2 classes of large genomic zones with a robust replication timing across cell lines that do not behave in the same way. The rst one corresponds to the large domains of early replication associated with a conserved RFD prole with a high density of replication initiation zones. The second class corresponds to robustly late replication domains with a variable RFD prole between cell line and a low density of initiation zones. Along the latter class, we visualized late replicating cell line specic initiation zones. Experiments demonstrated that the activity of replication origins depends on the chromatin environment[START_REF] Macalpine | Chromatin and DNA replication[END_REF].Late replicating regions are associated with a low density of transcribed genes and there is a strong correlation between chromatin modications associated with active transcription and early replication. So, we can hypothesize that these late initiation zones are associated with replication specic epigenetics marks, allowing the identication of the specic epigenetics mark related to replication initiation independently of transcription.

Figure A. 1 :Figure A. 2 :Figure A. 3 :

 123 Figure A.1: RFD proles correlation matrix between 2 technical replicates of IMR90 and 1 replicate of hTERT immoratilized IMR90 cell lines at 10 kb.

Figure A. 4 :

 4 Figure A.4: Detection of Stable and Variable regions. RFD correlation matrices of the 5Mb non overlapping window of chromosome 1. Each square represent a 5 Mb RFD correlation matrix among the 12 cell lines at 10 kb. Pearson correlation coecient values are color-coded from blue (0.4) to red (1) as in Figure 3.7. Order is upward then rightward from bottom left corner.

Figure A. 5 :

 5 Figure A.5: Overview of human chromosomes and their MCR levels. We map the chromosomes between 1 and 11 using the MCR at 100 kb. Colors code correspond to the 4 MCR levelas in Figure 3.18.

Figure A. 6 :

 6 Figure A.6: Overview of human chromosomes and their MCR levels. We map the chromosomes between 12 and 22 using the MCR at 100 kb. Colors code correspond to the 4 MCR levelas in Figure 3.18.

Figure A. 7 :

 7 Figure A.7: Number of common IZ with slope>2%RFD per kb at 30 kb resolution for all pairs of cell lines.The upper triangular part of the matrix represents the observed count of common origin with slope>2%RFD per kb between pair of cell lines in 30 kb windows. The lower triangular part of the matrix represents the expected count of common IZ computed >2%RFD per kb in each cell lines. Main diagonal is the number of detected IZ of slope>2%RFD per kb in each cell lines.

Figure A. 8 :Figure A. 9 :

 89 Figure A.8: Matrix of observed/expected numbers of common for each cell line pair for a threshold of 2%RFD per kb. Each square of this matrix represents the normalized count of common origin between pair of cell lines. Normalization was done by dividing the observed count of common IZ (upper part of the matrix in Figure A.7) by the expected one (lower part of the matrix in Figure A.7).

Figure A. 10 :

 10 Figure A.10: Matrix of observed/expected numbers of common for each cell line pair for a threshold of 1%RFD per kb. Each square of this matrix represents the normalized count of common origin between pair of cell lines. Normalization was done by dividing the observed count of common IZ (upper part of the matrix in Figure A.9) by the expected one (lower part of the matrix in Figure A.9).

Figure B. 1 :

 1 Figure B.1: Chromosomes 1 to 14 are shown, see Figure B.2.

Figure B. 2 :Figure B. 4 :Figure B. 5 :

 245 Figure B.2: Cell line classication based on correlations between replication and gene expression proles for each chromosome. Correlation matrices between RFD proles (CRFD ; A), RNA-seq (CFPKM ; B) and MRT proles (CMRT; C); Pearson correlation coecient values are colour-coded from blue (0.4) to red (1) using the colour bar on the right. A corresponding dendrogram representation of the hierarchical classication of cell lines is shown on top of each correlation matrix; ordinate is the correlation distance. Chromosomes 15 to 22 are shown.

Figure B. 7 :

 7 Figure B.7: Principal component of the gene expression proles for all RNA-seq sqmples as computed by DESeq2.

Figure C. 1 :Figure D. 1 :Figure D. 2 :

 112 Figure C.1: Molecular mechanism of the CML. A break occurs in the abl gene on chromosome 9 and another break in the bcr gene on chromosome 22. The abl gene is a proto-oncogene. The translocation causes the fusion of these two genes bcr/abl. By its fusion with another gene, the proto-oncogene abl is then transformed into an oncogene. Therefore a new mRNA is produced: bcr/abl-mRNA. A bcr/abl fusion protein is thus synthesized. It induces a tyrosine kinase activity above normal which stimulates an over proliferation of precursor bone marrow cells. (Adapted from http://www.embryology.ch/francais/kchromaber/popupchromaber/ 02abweichende/mfphilainterak/01.html)

Figure D. 6 :

 6 Figure D.6: Procedure to obtain P493 cell lines. A normal B lymphocyte was rst transfected with EBNA under the control of a promoter that is turned on in the presence of estrogen. It was further transfected by Myc under the control of a promoter that is turned o in presence of tetracycline. Cell line was grown in 3 conditions, with only EBNA expression for long term, with only Myc over-expression for long term or with both Myc over-expression and EBNA expression for long term.

Table 3 .

 3 1: Ok-seq sequencing statistics. Cell lines are listed according to their order in most gures in this thesis. Number of reads in each sample are in the second column, biological replicate number and sample name in column 3 and 4.

Table 3 .

 3 2:The number of detected origin in each cell lines. First column represent the cell lines. second one represent the count af all detected origin per cell lines. The third (resp. forth and fth) represent the count of origin according to the maximum slope 'MS'>0.5%RFD per kb (resp. >1%RFD per kb and > 2%RFD per kb).

maximum and local minimum of the RFD slope proles.

  .2).

	IB118	8.3	7.37	8.69		10.56	9.86		10.63	11.29
	TLSE19	7.03	5.96	6.87		7.65		7.49		8.06		7.92	6.88
	IMR90	7.4	5.84	6.71		7.26		7.02		7.6		7.55	6.76
	Hela	4.98	4.06	4.48		4.5		4.65		4.93		4.7	4.76
	K562	7.35	7.75	8.66		8.8		9.19		10.09	10.86	4.76	6.76	6.88	8.82
	TF1_BcrAbl_6M	7.79	7.46	8.79		9.59		12.88	14.02		10.86	4.7	7.55	7.92	11.29
	TF1_BcrAbl_1M	8.01	7.35	8.65		9.44		13.77			14.02	10.09	4.93	7.6	8.06	10.63
	TF1_GFP	7.47	6.95	8.28		9.06				13.77	12.88	9.19	4.65	7.02	7.49	9.86
	IARC385	9.3	9.98	11.84			9.06		9.44		9.59	8.8	4.5	7.26	7.65	10.56
	BL79	9.3	10.98		11.84	8.28		8.65		8.79	8.66	4.48	6.71	6.87	8.69
	Raji	8.07		10.98	9.98		6.95		7.35		7.46	7.75	4.06	5.84	5.96	7.37
	GM06990		8.07	9.3		9.3		7.47		8.01		7.79	7.35	4.98	7.4	7.03	8.3
	GM06990	Raji	BL79	IARC385	TF1_GFP	TF1_BcrAbl_1M	TF1_BcrAbl_6M	K562	Hela	IMR90	TLSE19	IB118

  's lymphoma cell line Raji (eciency loss > 0.5%RFD per kb) in relation to gene expression data in these cell lines allowed us to further quantify gene expression changes in regions of IZ eciency weakening. 119 IZ (26.3%) are not associated to transcription in either cell line (closest expressed protein coding gene at least 50kb away from the IZ). For the 333 remaining IZ, the closest expressed gene changes expression between GM06990 and Raji in 162 cases (|FPKM GM 06990 -FPKM Raji |>2), correspond to 116 gene repressions (FPKM GM 06990 > FPKM Raji ) and 46 gene activations. In other words, if 71.6% (116/162) of gene expression changes associated to IZ eciency losses are gene repression as expected for a co-regulation of replication and transcription, 74.3% ((452-116)/452) of the IZ eciency losses are not associated to transcriptional repression. Reproducing the analysis between closely cell lines TF1_GFP and TF1_BcrAbl_1M provided similar results with 86.0% (37/43) of IZ eciency losses not being associated to transcriptional repression, in particular because in this case gene expression changes are very few.

  progression system. This cancer is characterized by the presence of the Philadelphia chromosome that contains a fusion gene called BCR-ABL1 resulting from the combination of the BCR and ABL genes as shown in the FigureC.1[START_REF] Lozzio | Human chronic myelogenous leukemia cell-line with positive philadelphia chromosome[END_REF][START_REF] Kurzrock | The molecular genetics of philadelphia chromosomepositive leukemias[END_REF]. The BCR gene is located on the chromosome 22 and the

		A			B
		1.0
		0.8
	cdf	0.0 0.2 0.4 0.6		| RFD 1M GFP | 1% Largest 5% Largest 10% Largest 20% Largest All
		0.3	0.4	0.5	0.6
		GC content (200 kb)

Table 5 .

 5 ). Step 1 changes, whatever their direction, were most frequently 1: Database of replication initiation zone change of eciency. Summary of the manual annotation of RFD proles for changes in IZ eciency between TF1_GFP and TF1_BrcAbl_1M (top), TF1_BrcAbl_1M and TF1_BrcAbl_6M (middle) and, TF1_BrcAbl_6M and K562 cell lines (bottom). Note that for the latter case, only the 1027 loci presenting an IZ eciency change in at least one of the two rst comparisons were analysed. IZ change types (New, Enhanced, Weakened and Silenced) were organised to highlight the active or inactive status of IZ in each cell line. Figure 5.6: Persistence of replication initiation zones eciency changes in the CML progression model.(A) Distribution of IZ eciency changes (New, Enhanced, Weakened and Silenced) in Step 1 (black; n= 476) and Step 2 (yellow; n=774) of CML progression). (B) Evolution of Step 1 changes during Step 2, reported separately for each type of change. (C) Evolution of Step 2 changes in K562, reported as in (B). Error bars represent a 95% condence interval assuming data counts follow a Poisson distribution.

	Chapter V

Step 2 (n=88) were most frequently enhanced (n=50). In contrast, weakened IZs at Step 1 that changed again at Step 2 (n=131) were most frequently further weakened (n=89) or silenced (n=33). Therefore, there was a signicant tendency for IZ eciency changes during Step 2 to follow the same direction as observed during Step 1. These results demon-strate that BCR_ABL1 expression gradually alters the eciency of specic IZs in TF1 cells during long-term growth, resulting in a progressive destabilization of their replication program.

Table D .

 D 1: Description of the P493 cell line datasets.

	Short name	Expression of genes	Culture conditions
	P493_mEl	With only EBNA expression for long term (≥ 30 days)	+estrogen +tetracycline for ≥ 30 days
	P493_MEl	With both Myc over-expression and EBNA expression for long term (≥ 30	+estrogen -tetracycline for ≥ 30 days
		days)
	P493_Mel	With only Myc over-expression for long term (≥ 30 days)	-estrogen -tetracycline for ≥ 30 days
	P493_mEl_MEs With only EBNA expression for long	+estrogen +tetracycline for
		term (≥ 30 days ), then with both Myc over-expression and EBNA expression for short term (1 day)	≥ 30 days, then +estrogen -tetracycline for 1 day
	P493_mEl_Mes With only EBNA expression for long	+estrogen +tetracycline for
	term (≥ 30 days), then with only Myc over-expression for short term (1 day) P493_Mel_mEs With only Myc over-expression for long term (≥ 30 days), then with only EBNA	≥ 30 days, then -estrogen -tetracycline for 1 day -estrogen -tetracycline for ≥ 30 days, then +estrogen +tetracy-
		expression for short term (1 day)	cline for 1 day
	P493_MEl_mEs With both Myc over-expression and EBNA expression for long term (≥ 30	+estrogen -tetracycline for ≥ 30 days, then +estrogen +tetracy-
		days), then with only EBNA expression	cline for 1 day
		for short term (1 day)

The eukaryotes gather four kingdoms of the living world: animals, mushrooms, plants and protists which are characterized by the presence of a nucleus

The nucleosome is the fundamental unit of chromatin, consisting of a segment of DNA wound in sequence around eight histone protein cores.

CpG dinucleotide is a two-nucleotide DNA segment whose nucleic base sequence is CG.

A G-quadruplex (G4) is a four-strand secondary structure that can be adopted by nucleic acids rich in guanine.

15 25 35 45 55

75 85 95 105 115

Thus in summary, on the rst hand, we detected a coupling between transcriptional changes at some loci but initiation zones changes and on the second hand, no IZ eciency changes at other loci of dierentially expression.

This protein is able to induce transcription through a pyrimidine-rich promoter element.

20170317_P493_Mel R20171011_P493_Mel_S5 R20171206_P493_Mel_S3 R20171206_P493_Mel_S5 R20171011_P493_mEl_Mes_S6 R20171011_P493_mEl_S1 R20171011_P493_mEl_S2 R20171206_P493_mEl_S1 20160421_P493_MEl_mEs R20171011_P493_MEl_S4 R20171206_P493_MEl_S2 R20171206_P493_MEl_S4 R20171206_P493_MEl_S6 20170317_P493_MEl GM06990 20170317_P493_Mel R20171011_P493_Mel_S5 R20171206_P493_Mel_S3 R20171206_P493_Mel_S5 R20171011_P493_mEl_Mes_S6 R20171011_P493_mEl_S1 R20171011_P493_mEl_S2 R20171206_P493_mEl_S1 20160421_P493_MEl_mEs R20171011_P493_MEl_S4 R20171206_P493_MEl_S2 R20171206_P493_MEl_S4 R20171206_P493_MEl_S6 20170317_P493_MEl

Supplementary gures for Chapter IV