Temporal association of ORCA/LRWD1 to late-firing origins during G1 dictates heterochromatin replication and organization

Chromatin immunoprecipitation (ChIP) and ChIP-seq

Detailed ChIP protocol is described in supplementary material and methods.

Construction of ChIP-Seq libraries and sequencing on the HiSeq2500 was carried out at the Roy J. Carver Biotechnology Center, University of Illinois at Urbana-Champaign (UIUC). The libraries were constructed with the Kapa Hyper Prep Kit from Kapa Biosystems (Kapa Biosystems, MA, USA). Briefly, DNAs (10 ng) were blunt-ended, 3΄-end A-tailed and ligated to indexed adaptors with 6nt barcodes. Adaptored DNA were amplified by PCR to selectively enrich for those fragments that have adapters on both ends. Amplification was carried out for 10 cycles with the Kapa HiFi polymerase (Kapa Biosystems, Woburn, MA, USA) to reduce the likeliness of multiple identical reads due to preferential amplification. The final libraries were run on Agilent bioanalyzer DNA high-sensitivity chip (Agilent, Santa Clara, CA, USA) to determine the average fragment size and to confirm the presence of DNA of the expected size range. They were also quantitated by qPCR on a BioRad CFX Connect Real-Time System (Bio-Rad Laboratories, Inc. CA, USA) prior to pooling and sequencing.

The ChIP-Seq libraries were pooled in equimolar concentration based on the qPCR concentration and sequenced on four lanes on an Illumina HiSeq2500 using HiSeq SBS sequencing kits version 4. The raw .bcl files were converted into demultiplexed compressed fastq files using the bcl2fastq 1.8.2 conversion software (Illumina). Independently, ChIP-Seq for the 1.5 and 3 h time points was performed in duplicates with two independent biological replicates using the HiSeq Rapid choice method.

Nascent strand sequencing

Replication origins were mapped using nascent strand sequencing (19,35). In short, the nascent strand samples were prepared by collecting λ exonuclease resistant, 0.5–2.5 kb fragments of single stranded DNA from asynchronous U2OS cell populations. High molecular weight genomic DNA was prepared in parallel from the same population of cells, sonicated and sequenced. All samples were subject to massively parallel sequencing using standard Illumina protocols (36).

Following sequencing, Genomatix (https://www.genomatix.de/) was used to call peaks against appropriate genomic DNA or input controls. Either the MACS algorithm or the SICER algorithm was used. SICER was used to call peaks for the nascent strand sequencing with the following parameters: redundancy threshold = 2, window size = 200, fragment size = 150, gap size = 600, FDR = 0.01, P-value = 200. MACS was used to call peaks for the ORCA ChIP-Seq with the following parameters: Tag size = 20, q-value = 0.01, bandwidth = 300, lower and upper limit to mold for (model) = 10–30, redundancy threshold = auto. Bed file peaks identified by SICER or MACS were visualized along with the tiled original data using the Integrative Genome Viewer (https://www.broadinstitute.org/igv/).

BED file comparisons

A custom script (available upon request) was used to create a list of sequences from a reference file (input 1) that are or are not found within a user-defined distance of sequences in the comparator file (input 2). The BedIntersect script outputs the peaks from the reference file that overlap within 2 kb of the comparator, while the BedSubtract script outputs the peaks that are found in the reference and not in the comparator. Series of intersections and subtractions were performed to identify the extent of colocalization among ORCA, methylated CpGs, histone modifications and replication origins. Whole-genome co-localization was visualized and quantified using ColoWeb (37) and verified with Genomatix.

Methylated DNA immunoprecipitation (MeDIP) and MeDIP-seq

Detailed MeDIP protocol is described in supplementary material and methods.

For MeDIP-seq, sonicated genomic DNA (6 μg) was used to construct libraries using the Kapa Hyper Prep Kit (Kapa Biosystems, MA) prior to performing MeDIP. Briefly, 6 μg of sonicated DNA were blunt-ended, 3΄-end A-tailed and ligated to indexed adaptors with 6nt barcodes. Adaptored DNA was then subjected to MeDIP as described above. Five nanogram of pull down DNA were amplified by PCR to selectively enrich for those fragments that have adapters on both ends. Amplification was carried out for 12 cycles with the Kapa HiFi polymerase (Kapa Biosystems, Woburn, MA, USA) to reduce the likeliness of multiple identical reads due to preferential amplification. The final libraries were run on Agilent bioanalyzer DNA high-sensitivity chip (Agilent, Santa Clara, CA, USA) to determine the average fragment size and to confirm the presence of DNA of the expected size range. They were also quantitated by qPCR on a BioRad CFX Connect Real-Time System (Bio-Rad Laboratories, Inc., CA, USA) prior to pooling and sequencing.

The MeDIP libraries were pooled in equimolar concentration based on the qPCR concentration and sequenced on two lanes on an Illumina HiSeq2500 using HiSeq SBS sequencing kits version 4. The raw .bcl files were converted into demultiplexed compressed fastq files using the bcl2fastq 1.8.2 conversion software (Illumina).

MeDIP-seq data analysis

The MeDIP-seq sequencing reads were aligned to the GRCh37/hg19_noMask version of the Human genome using Bowtie 0.12.7 (38). The single-end MeDIP-seq sequencing reads were aligned with the following Bowtie parameters: -n 2 -l 20 -M 1 –tryhard. MeDIP-seq peaks were called using Model-based Analysis of ChIP-Seq (MACS) with default parameters (39). Input sequencing libraries were used to control for sequencing specific biases.

Data access

All sequencing data have been submitted to NCBI Gene Expression Omnibus (GEO) (http://www.ncbi.nlm.nih.gov/geo/). The accession numbers are: H3K9me3 ChIP-seq (GSE68129); ORCA ChIP-seq and MeDIP-seq in U2OS cells (GSE81165); U2OS nascent strand-seq (GSE80391).

Binding of ORCA on chromatin is dynamically regulated during G1

To understand how ORCA regulates heterochromatin replication and/or organization, we evaluated the genome-wide distribution of ORCA using chromatin immunoprecipitation (ChIP) in U2OS cell line stably expressing HA-ORCA (Supplementary Figure S1A). We performed ChIP in cells synchronized at G1 phase, a stage when ORC loading and the licensing process are established (Figure 1). We found that ORCA showed dynamic binding pattern at several genomic loci at different stages of G1. Previous work has shown that the largest ORC subunit, Orc1 protein shows spatio-temporal dynamics during G1 (40). Since the protein levels of ORCA are highest during G1, we conducted ORCA ChIP at different points in G1. Cells were synchronized at prometaphase by nocodazole, and the mitotic cells were collected by shake-off. The cells were released from the arrest and collected for ChIP at various time points within G1 [1.5 h (early G1), 3 h (mid G1) and 5 h (late G1)]. The arrest during mitosis and the synchronous release was monitored by flow cytometry, immunoblot analyses of G1 and S phase markers (Supplementary Figure S1B) and immunofluorescence to evaluate BrdU incorporation (not shown). ChIP was performed using HA antibody followed by deep sequencing at early, mid, and late G1 phase. Our ChIP-seq result showed that both the binding pattern and the number of ORCA binding sites changed as cells progress through G1 phase. We obtained ∼10 000 peaks at early and mid G1 and only 1347 peaks at late G1 (Figure 1A). This change in the number of ORCA binding sites during G1 is consistent with the fact that ORCA protein level peaks at G1 and drops at the G1/S boundary (15). We also found that about 55% of the peaks at early and mid G1 are not found at late G1, whereas most of the peaks at late G1 are also present in the previous two time points (Figure 1A and B). A subset of regions is bound by ORCA throughout G1, which contribute to most of the ORCA binding regions we have identified at the late G1 time point (Figure 1A). Representative snapshots of regions with different ORCA binding pattern throughout G1 are shown (Figure 1Ba, Bc, Be). Binding pattern of ORCA for c-Fos locus as a negative control is shown in Supplementary Figure S1C. These results indicate that the binding of ORCA onto chromatin during G1 phase is a temporally regulated process and is dynamic. The temporal binding of ORCA at specific regions was also validated by ChIP-qPCR (Figure 1Bb, Bd, Bf). Similar temporal ORCA binding on chromatin was also observed by ChIP-qPCR in wild-type U2OS cells using ORCA antibody (Supplementary Figure S1Da and Db).

ORCA binding regions contained ∼32–34% intergenic regions, with the other areas either genic or partial (partly contained in genic regions) (Supplementary Figure S1E). Since the ChIP-Seq analysis is limited to single copy regions and does not include repetitive loci, which comprise a large fraction (>50%) of the human genome, the actual abundance of ORCA binding sites intergenic regions is likely to be higher. Remarkably however, ORCA peaks are significantly enriched in promoter regions (Supplementary Figure S1E), consistent with their colocalization with replication origins (see below).

ORCA binding regions strongly co-localize with replication origins

We have previously shown that ORCA interacts with the ORC complex and stabilizes ORC on the chromatin (15). Further, loss of ORCA in primary fibroblasts resulted in cells arresting in G1 with reduced MCM loading, consistent with a role for ORCA in replication origin licensing. In addition to its role in replication initiation, ORCA also regulates heterochromatin organization in a way that is independent of its role in replication initiation (33). To address how ORCA might regulate replication initiation, we asked if the ORCA binding sites during G1 represent replication origins. We mapped the replication initiation sites in U2OS cells by nascent strand abundance analysis (35). Our analysis showed that there is strong co-localization between ORCA-binding sites and replication initiation sites. At all three time points tested, most of the ORCA binding sites (65–85%) were found near (within 2 kb) replication origins (Figure 2A and B). However, these ORCA binding sites represent only 11.5% of the total mapped origins. ORCA binding sites, however, did not show strong co-localization with a random generated genome sequence of the same size distribution and GC content as the origin sequences (Figure 2B), suggesting that the majority of ORCA binding sites identified in our ChIP-seq analysis indeed represent replication origins. We observed that at each time point, ORCA specifies a different subset of origins, indicating that ORCA associates with distinct origins during G1 in a dynamic manner. In fact, at early and mid G1, more than 50% of the origins bound by ORCA is not found at late G1, while a subset of origins are bound by ORCA throughout G1 (Figure 2C and D). In late G1 phase, most of the origins bound by ORCA are present in the two earlier time points and very few origins are bound by ORCA only at late G1 (Figure 2C). These results indicate that the specification of replication origins by ORCA/ORC during G1 may be more dynamic than previously thought. Once an origin is specified by ORCA, it may enable the licensing process and once the licensing process is finished, ORCA may dissociate from its initial loading site.

ORCA binding regions are enriched for H3K9me3 and methyl-CpG marks

We previously demonstrated that ORCA localizes at repressive chromatin structures and regulates heterochromatin organization (15,33). We mapped the genome-wide distribution of repressive marks including H3K9me3 and DNA methylation in human U2OS cells. By aligning the data obtained from the ORCA ChIP-seq with that of H3K9me3 ChIP-seq, we found that there is strong co-localization between ORCA-binding regions and H3K9me3-containing regions (Figure 3Aa and c, Supplementary Figure S2A). Interestingly, the percentage of ORCA peaks that co-localize with H3K9me3 increases slightly as cells progress through G1 (from 45% to 65%, Figure 3Ab) consistent with our cell biological data that ORCA levels decrease at the end of G1 but whatever remains localizes to heterochromatic structures.

To address ORCA's localization in relation to DNA methylation, we carried out methylated DNA immunoprecipitation (MeDIP) with an antibody that specifically recognizes methyl-cytosine. We performed MeDIP followed by deep sequencing in U2OS cells to determine the genome wide localization of methyl-CpG sites. We aligned the MeDIP-seq result with our ORCA ChIP-seq results. Similar to H3K9me3, we also found strong co-localization between ORCA-binding regions and methyl-CpG sites, and the percentage of ORCA peaks that co-localize with methyl-CpG sites also increases slightly as cells progress through G1 (from 49% to 61%, Figure 3Ba and Bb).

ChIP-Seq was performed in two biological replicates for the 1.5 and 3 h time points. In both sets of biological replicates we observed high concordance between the duplicates collected at the same time but low concordance between the two time points. Both biological replicates demonstrated a highly significant colocalization with replication origins and H3K9me3 when compared to randomized controls. These data support our conclusion that ORCA associates with replication origins and H3K9 methylated regions and that ORCA binding undergoes temporal relocation during G1.

ORCA regulates H3K9me3 and methyl-CpG marks at its binding sites

We have previously shown that the loss of ORCA causes dramatic reduction of H3K9me3 at 18% of the H3K9me3 peaks (33). However, we were unable to pinpoint why only specific H3K9me3 sites were altered. To address if ORCA regulates H3K9me3 at its binding sites, we aligned the H3K9me3-containing ORCA peaks with the 18% H3K9me3 peaks (ORCA-dependent H3K9me3 peaks) that showed dramatic reduction upon ORCA depletion. Our results showed that most of the H3K9me3-containing ORCA-binding sites (∼80%) showed H3K9me3 reduction upon ORCA depletion (Figure 3C, Supplementary Figure S2B), indicating that ORCA binds to H3K9me3-containing regions and regulates the local repressive marks. Two snapshots showing the loss of H3K9me3 upon ORCA-depletion occurring at the regions showing co-localization of ORCA and H3K9me3 are provided (Supplementary Figure S2B). Adjacent regions enriched for H3K9me3 (without ORCA binding) remain unaffected upon the loss of ORCA (Supplementary Figure S2B). However, what dictates ORCA-binding to certain repressive sites and not others remains to be determined. Future studies will determine the potential involvement of other factors, chromatin modifications, or DNA methylation in the differential binding of ORCA to repressive marks.

Previous work demonstrated that ORCA binds to nucleosomes most efficiently in the presence of methylated DNA (31). We therefore evaluated if ORCA also regulates DNA methylation status at its binding sites. We carried out MeDIP followed by qPCR in control and ORCA-depleted cells. Our MeDIP-qPCR results showed a significant decrease in DNA methylation at several ORCA-binding sites (the ones that are also enriched for H3K9me3) upon ORCA depletion (Figure 3D). Interestingly, we observed about a 50% reduction in DNA methylation level near the transcription start sites of ribosomal RNA (Figure 3D). The DNA methylation status at the rDNA locus is known to control rRNA transcription (41). We therefore determined if the methylation change caused by ORCA depletion changed the rRNA transcription status. Our transcription analysis by qPCR showed an increase in 45S pre-rRNA level in ORCA-depleted cells (Supplementary Figure S2C), indicating that the loss of DNA methylation at rDNA locus by ORCA depletion also causes a change in rRNA transcription. Previously, it has been shown in mouse cells that loss of ORCA leads to increased major satellite repeat transcription (42). Here, it appears that ORCA also regulates transcription from certain repeat regions in human cells.

ORCA-bound origins are enriched for H3K9me3 and methyl-CpG marks

Our analysis showed that most of the ORCA-binding regions co-localize with replication origins. However, there appear to be a large number of origins that are not bound by ORCA. This is similar to what has been reported previously for Orc1 binding to human genome (21). We therefore determined the chromatin feature of a subset of distinct origins that are occupied by ORCA. Replication origins are generally enriched at euchromatic regions (12). Further, our analysis also showed that replication initiation events do not show very strong co-localization with H3K9me3 and methylated CpG sites (Figure 4A). In fact, among the total replication initiation sites, only a small subset is found near H3K9me3 (20%) and 40% of the origins are found near methyl-CpG site (Figure 4A). However, ORCA-occupied origins showed strong co-localization with H3K9me3 and methyl-CpG sites (Figure 4B and C). At all the three time points in G1, a large portion of the ORCA-occupied origins are found near H3K9me3 and methyl-CpG sites (Figure 4B and C). A snapshot of an ORCA-bound origin enriched for repressive environment is shown (Figure 4D). These results suggest that ORCA may specify a subset of replication origins with repressive chromatin marks, and facilitate the loading of other pre-RC components at these origins.

Our observation that ORCA-bound regions is enriched for repressive chromatin marks indicates that these regions may be late-replicating regions. Using the recently published U2OS replication timing profile (43), we stratified the genome into regions that replicate during early, early-mid, mid-late and late S phase. Our analysis showed that regions that are bound by ORCA replicated predominantly during late S phase (Supplementary Figure S2D). The numbers of binding regions that intersect with late replicating regions were higher during early G1 (1.5 and 3 hours into G1) but the enrichment was also evident at the 5 hours time point, albeit with fewer peaks overall (73.0% of ORCA 1.5 h; 70.7% of ORCA 3 h and 67.1% of ORCA 5 h, found in late replicating regions). This observation is in line with the notion that repressive regions generally replicate late during S phase. Our analysis showed that H3K9me3 containing ORCA peaks primarily colocalize with origins (Supplementary Figure S2E). This is consistent with the fact that ORCA-bound origins are late-firing origins.

ORCA interacts with methylated DNA sequence in vitro and associates with DNA methyltransferases

Our sequencing analyses showed that a large number of ORCA peaks co-localize with H3K9me3 and methyl-CpG peak in vivo. Others and we have previously shown that ORCA could interact with both H3K9 methyltransferases and methylated H3K9 peptide in vitro (33,42). We now addressed if ORCA can bind methylated DNA as well. To address this, we performed electrophoretic mobility shift assay (EMSA) using bacterially purified His-SUMO-ORCA protein. To test whether ORCA interacts with methylated DNA directly, we used a double-stranded DNA probe with 12 CpG sites, with the cytosine being either un-modified or methylated. When we incubated the methylated probe with His-SUMO-ORCA but not His-SUMO, we observed a dramatic shift in mobility, indicating a strong interaction between ORCA and the methylated DNA probe (Figure 5A). Incubation of an un-methylated probe with ORCA, however, did not show an obvious shift, indicating ORCA specifically recognizes and interacts with methylated DNA. We then asked if ORCA only recognizes methyl-cytosine when the CpG sites on both strands are methylated. To address this, we generated DNA substrates that are unmethylated, hemi-methylated (one strand) or methylated (both the strands) and performed EMSA experiment with purified ORCA protein. The result showed that ORCA interacts with both the hemi-methylated and fully methylated probe but not with the un-methylated probe (Figure 5B), suggesting that ORCA may be recognizing the methyl-group on methyl-cytosine. We also performed a similar experiment with another pre-RC protein, Orc1, and did not observe any binding to methylated or unmethylated DNA substrate (Supplementary Figure S3A).

We found that ORCA binds to repressive histone marks and also associates with the enzymes (HKMTs) that catalyze these marks (33). Since we found that ORCA preferentially binds to methylated sequences, we examined if ORCA also interacts with the DNA methylation machinery, DNA methyltransferases (DNMTs). Mammalian cells have three active DNA methyltransferases: DNMT1, DNMT3a and DNMT3b (44). To address whether ORCA interacts with DNA methyltransferases, we performed co-immunoprecipitation (co-IP) experiments. By performing co-IP using T7 antibody in U2OS cells expressing Myc-DNMT1 with or without T7-ORCA, we found that ORCA interacts with DNMT1 (Figure 5C). In a U2OS cell line stably expressing HA-ORCA, we performed IP using HA antibody and the result revealed robust interaction between HA-ORCA and endogenous DNMT3a (Figure 5D). Similarly, endogenous ORCA was found to interact with Myc-DNMT3b (Figure 5E). Reciprocal co-IPs confirmed these interactions (Supplementary Figure S3B–D). To rule out the possibility that DNA or RNA mediated these interactions, we have also performed co-IP experiments in the presence of Ethidium Bromide (EtBr) or RNase A. Interactions between DNMT1, DNMT3a, DNMTb and ORCA were still detected under these conditions (Supplementary Figure S3E–I). These results demonstrate that ORCA associates with methylated DNA and also with the enzymatic machinery catalyzing these marks.

Repressive marks are required for ORCA binding to chromatin

We have observed that ORCA binding regions are enriched for repressive marks including H3K9me3 and methylated CpG sequence in vivo. We observe that ORCA directly interacts with both methylated H3K9 peptides and methylated DNA sequences in vitro. We reasoned that these repressive marks could also be influencing ORCA recruitment onto chromatin. To test this hypothesis, we depleted H3K9 methyltransferases SUV39H1 and SUV39H2 by siRNA (Supplementary Figure S3J). SUV39H1 and H2 are responsible for establishing H3K9me3 and depletion of these two enzymes leads to a significant loss of H3K9me3 mark (Supplementary Figure S3J and S3K). It has been shown that H3K9me3 and DNA methylation are dependent on each other (45), and we found that depletion of SUV39H1 and H2 also led to a reduction of DNA methylation level at several loci (Supplementary Figure S3L). We then carried out ORCA ChIP in SUV39-depleted cells and found that ORCA binding at several genomic loci is significantly reduced upon SUV39 depletion (Figure 5F), indicating that these repressive chromatin marks are required for ORCA binding at these regions. The total protein level of ORCA did not change significantly in SUV39-depleted cells (Supplementary Figure S3J). Thus, ORCA could be recruited to its binding sites by directly recognizing repressive chromatin marks and hence recruit other factors to its binding sites for its replication and/or heterochromatin functions.