The effect of Nipped-B-like (Nipbl) haploinsufficiency on genome-wide cohesin binding and target gene expression: modeling Cornelia de Lange syndrome

Cells and antibodies

Mouse embryonic fibroblasts (MEFs) derived from E15.5 wild type and Nipbl mutant embryos were used as described previously [19]. In brief, mice heterozygous for Nipbl mutation were generated (Nipbl+/−) from gene-trap-inserted ES cells. This mutation resulted in a net 30–50% decrease in Nipbl transcripts in the mice, along with many phenotype characteristics of human CdLS patients [19]. Wild type and mutant MEF cell lines derived from the siblings were cultured at 37 °C and 5% CO₂ in DMEM (Gibco) supplemented with 10% fetal bovine serum and penicillin-streptomycin (50 U/mL). Antibodies specific for hSMC1 and Rad21 were previously described [46]. Rabbit polyclonal antibody specific for the NIPBL protein was raised against a bacterially expressed recombinant polypeptide corresponding to the C-terminal fragment of NIPBL isoform A (NP_597677.2) (amino acids 2429–2804) [45]. Anti-histone H3 rabbit polyclonal antibody was from Abcam (ab1791).

ChIP-sequencing (ChIP-seq) and ChIP-PCR

Peak finding

Peaks were called using AREM (Aligning ChIP-seq Reads using Expectation Maximization) as previously described [50]. AREM incorporates sequences with one or many mappings to call peaks as opposed to using only uniquely mapping reads, allowing one to call peaks normally missed due to repetitive sequence. Since many peaks for Rad21 as well as CTCF can be found in repetitive sequence [50, 51], we used a mixture model to describe the data, assuming K + 1 clusters of sequences (K peaks and background). Maximum likelihood is used to estimate the locations of enrichment, with the read alignment probabilities iteratively updated using EM. Final peaks are called for each window assuming a Poisson distribution, calculating a p value for each sequence cluster. The false discovery rate for all peaks was determined relative to the pre-immune sample, with EM performed independently for the pre-immune sample as well. Full algorithm details are available, including a systematic comparison to other common peak callers such as SICER and MACS [50]. Overlap between peaks and genomic regions of interest were generated using Perl and Python scripts as well as pybedtools [52, 53]. Figures were generated using the R statistical package [54]. Visualization of sequence pileup utilized the UCSC Genome Browser [55, 56].

Motif analysis

De novo motif discovery was performed using Multiple Expectation maximization for Motif Elicitation (MEME) version 6.1 [57]. Input sequences were limited to 200 bp in length surrounding the summit of any given peak, and the number reduced to 1000 randomly sampled sequences from the set of all peak sequences. Motif searches for known motifs were performed by calculation of a log-odds ratio contrasting the position weight matrix with the background nucleotide frequency. Baseline values were determined from calculations across randomly selected regions of the genome. Randomly selected 200-bp genomic regions were used to calculate a false discovery rate (FDR) at several position weight matrix (PWM) score thresholds. We chose the motif-calling score threshold corresponding to a 4.7% FDR. The p values were derived for the number of matches above the z-score threshold relative to the background using a hypergeometric test.

Expression data analysis

Affymetrix MOE430A 2.0 array data for mouse embryonic fibroblasts (10 data sets for the wild type and nine for Nipbl+/− mutant MEFS) were previously published [19]. Expression data were filtered for probe sets with values below 300 and above 20,000, with the remainder used for downstream analysis. Differential expression and associated p values were determined using Cyber-t, which uses a modified t test statistic [58]. Multiple hypothesis testing correction was performed using a permutation test with 1000 permutations of the sample data. Probe sets were collapsed into genes by taking the median value across all probe sets representing a particular gene. Raw expression values for each gene are represented as a z-score, which denotes the number of standard deviations that value is away from the mean value across all genes. Gene ontology analysis was performed using PANTHER [59, 60] with a cutoff of p < 0.05.

KS test

Genes were sorted by their fold-change, and any adjacent ChIP-binding sites were identified. We performed a Kolmogorov-Smirnov (KS) test comparing the expression-sorted ChIP binding presence vs. a uniform distribution of binding sites, similar to Gene Set Enrichment Analysis [61]. If ChIP binding significantly correlates with the gene expression fold-change, the KS statistic, d, will also have significant, non-zero magnitude. To better visualize the KS test, we plotted the difference between the presence of cohesin binding at (expression-sorted) genes in Fig. 5. The x axis of this figure is the (fold-change-based) gene rank, and the y axis is the KS statistic d, which behaves like a running enrichment score and is higher (lower) when binding sites co-occur more (less) often than expected if there were no correlation between ChIP binding and expression fold-change. The KS test uses only the d with the highest magnitude, which is indicated in the plots by a vertical red line. To better visualize ChIP binding presence, we further plot an x-mirrored density of peak presence at the top of each plot; the gray “beanplot” [62] at the top of the plots are larger when many of the genes have adjacent ChIP-binding sites.

siRNA depletion

Wild type MEFs were transfected using HiPerFect (Qiagen) following the manufacturer’s protocol with 10 mM small interfering RNA (siRNA). A mixture of 30 μl HiPerFect, 3 μl of 20 μM siRNA, and 150 μl DMEM was incubated for 10 min and added to 2 × 10⁶ cells in 4 ml DMEM. After 6 h, 4 ml fresh DMEM with 10% FBS was added. Transfection was repeated the next day. Cells were harvested 48 h after the first transfection. SiRNAs against Nipbl (Nipbl-1: 5′-GTGGTCGTTACCGAAACCGAA-3′; Nipbl-2: 5′-AAGGCAGTACTTAGACTTTAA-3′) and Rad21 (5′-CTCGAGAATGGTAATTGTATA-3′) were made by Qiagen. AllStars Negative Control siRNA was obtained from Qiagen.

RT-q-PCR

Total RNA was extracted using the Qiagen RNeasy Plus kit. First-strand cDNA synthesis was performed with SuperScript II (Invitrogen). Q-PCR was performed using the iCycler iQ Real-time PCR detection system (Bio-Rad) with iQ SYBR Green Supermix (Bio-Rad). Values were generated based on Ct and normalized to control gene Rnh1. PCR primers specific for major satellite, minor satellite, rDNA, and SINE B1 were previously described [48]. Other unique primers are listed in Table 1. The RT-qPCR analyses of the wild type and mutant cells were done with two biological replicates with consistent results. The gene expression changes after siRNA treatment were evaluated with two to three biological replicates with similar results.

3C analysis

The chromosome conformation capture (3C) protocol was performed as described [35]. Approximately 1 × 10⁷ cells were crosslinked with 1% formaldehyde at 37 °C for 10 min. Crosslinking was stopped by adding glycine to a final concentration of 0.125 M. Cells were centrifuged and lysed on ice for 10 min. Nuclei were washed with 500 μl of 1.2× restriction enzyme buffer and resuspended with another 500 μl of 1.2× restriction enzyme buffer with 0.3% SDS and incubated at 37 °C for 1 h. Triton X-100 was added to 2% and incubated for another 1 h. 800 U of restriction enzyme (HindIII New England Biolabs) was added and incubated overnight at 37 °C. The digestion was heat-inactivated the next day with 1.6% SDS at 65 °C for 25 min. The digested nuclei were added into a 7 ml 1× ligation buffer with 1% Triton X-100, followed by 1-h incubation at 37 °C. T4 DNA ligase (2000 U) (New England Biolabs) was added and incubated for 4 h at 16 °C followed by 30 min at room temperature. Proteinase K (300 μg) was added, and the sample was reverse-crosslinked at 65 °C overnight. Qiagen Gel Purification Kits were used to purify DNA. Approximately 250 ng of template was used for each PCR reaction. PCR products were run on 2% agarose gels with SYBRSafe (Invitrogen), visualized on a Fujifilm LAS-4000 imaging system and quantified using Multigauge (Fujifilm).

To calculate interaction frequencies, 3C products were normalized to the constitutive interaction at the excision repair cross-complementing rodent repair deficiency, complementation group 3 (ercc3) locus [63, 64], which is unaffected in mutant MEFs. A control template was made to control for primer efficiencies locus-wide as described [65]. PCR fragments spanning the restriction sites examined were gel purified, and equimolar amounts were mixed (roughly 15 μg total) and digested with 600 U restriction enzyme overnight and subsequently ligated at a high DNA concentration (> 300 ng/μl). The template was purified with the Qiagen PCR Purification Kit and mixed with an equal amount of digested and ligated genomic DNA. Two hundred fifty nanograms of the resulting control template was used for each PCR for normalization against PCR primer efficiencies. Two biological replicates with three technical replicates each were analyzed for both wild type and mutant cells and for control and Nipbl siRNA-treated cells, which yielded consistent results.

Nipbl haploinsufficiency leads to a global reduction of cohesin binding to its binding sites

In order to investigate how Nipbl haploinsufficiency leads to CdLS, cohesin binding was examined genome-wide by ChIP-seq analyses using antibody specific for the cohesin subunit Rad21, in wild type and Nipbl+/− mutant MEFs derived from E15.5 embryos [19] (Fig. 1a). MEFs derived from five wild type and five mutant pups from two litters were combined to obtain sufficient chromatin samples for ChIP-seq analysis. Nipbl+/− mutant MEFs express approximately 30–40% less Nipbl compared to wild type MEFs [19] (Table 2). MEFs from this embryonic stage were chosen in order to match with a previous expression microarray study, because they are relatively free of secondary effects caused by Nipbl mutation-induced developmental abnormalities compared to embryonic tissue [19]. Consistent with this, there is no noticeable difference in growth rate and cell morphology between normal and mutant MEFs [19]. This particular anti-Rad21 antibody was used previously for ChIP analysis and was shown to identify holo-cohesin complex binding sites [30, 35, 45, 66]. This is consistent with the close correlation of the presence of other cohesin subunits at identified Rad21-binding sites [67] (Fig. 1b).

Gene	Nipbl+/− mutant	Nipbl siRNA	Rad21 siRNA
Nipbl	0.68 ± 0.003	0.68 ± 0.001	1.04 ± 0.051
Rad21	0.94 ± 0.021	0.99 ± 0.021	0.26 ± 0.018
CTCF	0.95 ± 0.050	0.96 ± 0.066	0.84 ± 0.074

Cohesin-binding sites were identified using AREM [50], with a significance cutoff based on a p value less than 1 × 10⁻⁴, resulting in a FDR below 3.0% (Fig. 1a). Cohesin-binding peaks ranged from ~ 200 bp to ~ 6 kb in size with the majority less than 1 kb in both wild type and mutant cells (median value of 499 bp in wild type and 481 bp in mutant cells) (Fig. 1c). Approximately 35% fewer cohesin-binding sites were found in Nipbl+/− mutant MEFs compared to the wild type MEFs (Fig. 1a). This is not due to variability in sample preparation since no significant difference in the histone H3 ChIP-seq was observed between the wild type and mutant cell samples (R value = 0.96) (Fig. 1d). Since the total read number for mutant ChIP-seq was ~ 15% less than for wild type ChIP-seq (Fig. 1a), we examined whether the difference was in part due to a difference in the number of total read sequences between the two Rad21 ChIP samples. To address this, we randomly removed reads from the wild type sample to match the number of reads in the mutant sample and ran the peak discovery algorithm again on the reduced wild type read set. This was repeated 1000 times. We found that the wild type sample still yielded ~ 39% more peaks than the mutant, indicating that identification of more peaks in the wild type sample is not due to a difference in the numbers of total read sequences (Fig. 1e). Thus, cohesin appears to bind to fewer binding sites in Nipbl haploinsufficient cells.

The above results might suggest that a significant number of binding sites are unique to the wild type cells (Fig. 1a). When we compared the raw number of reads located within wild type peaks and the corresponding regions in mutant MEFs, however, we noted a reduced, rather than a complete absence of, cohesin binding in mutant cells (Fig. 1f). Those regions in mutant cells corresponding to the “WT only” regions consistently contain one to three tags in a given window, which are below the peak cutoff. However, the signals are significant compared to the negative control of preimmune IgG (Fig. 1f). Furthermore, even for those sites that are apparently common between the control and mutant MEFs, the binding signals appear to be weaker in mutant cells (Fig. 1f). To validate this observation, we segmented the genome into nonoverlapping 100 bp bins and plotted a histogram of the log ratios of read counts between the wild type and mutant samples in each bin, with read counts normalized using reads per kilobase per million total reads (RPKM) [68]. The plot indicates that the read counts for the mutant bins are generally less than those for the wild type bins, even for the binding sites common to both wild type and mutant cells (Fig. 1g). Signal intensity profiles of the Rad21 ChIP-seq in the selected gene regions also show a general decrease of Rad21 binding at its binding sites in Nipbl+/− MEFs compared to the control MEFs (see Fig. 6b). Decreased cohesin binding was further confirmed by manual ChIP-qPCR analysis of individual cohesin-binding sites using at least three independent control and mutant MEF samples supporting the reproducibility of the results (see Fig. 3). Decreased cohesin binding was also observed at additional specific genomic regions in Nipbl+/− MEFs [69]. Taken together, the results indicate that cohesin binding is generally decreased at its binding sites found in wild type MEFs, rather than re-distributed, in mutant MEFs.

The relationship of cohesin-binding sites with CTCF-binding sites and CTCF motifs

It has been reported that cohesin binding significantly overlaps with CTCF sites and depends on CTCF [30, 31]. A study in mouse embryonic stem cells (mESCs) showed, however, that there is only a limited overlap between CTCF- and Nipbl-bound cohesin sites, suggesting that there are two categories of cohesin-binding sites and the latter may be particularly important for gene activation [40]. Other studies also revealed that ~ 20–30% of cohesin sites in different human cancer cell lines and up to ~ 50% of cohesin sites in mouse liver appear to be CTCF-free [42, 43]. Some of these non-CTCF sites overlap with sequence-specific transcription factor binding sites in a cell type-specific manner, highlighting the apparent significance of CTCF-free cohesin sites in cell type-specific gene expression [42, 43]. De novo motif discovery by MEME identified the CTCF motif to be the only significant motif associated with cohesin-binding sites in our MEFs (Fig. 2a). Comparing our cohesin peaks with experimentally determined CTCF-binding peaks in MEFs [40], we found that approximately two thirds of cohesin-binding sites detected by Rad21 ChIP overlapped CTCF-binding sites (Fig. 2b). This is comparable with what was initially observed in mouse lymphocytes [30] and HeLa cells [31] using antibodies against multiple cohesin subunits. In contrast to recent studies reporting that almost all the CTCF-binding sites overlap with cohesin [43], our results show that less than 60% of CTCF-binding sites are co-occupied with cohesin (Fig. 2b). This is consistent with the fact that CTCF binds and functions independently of cohesin at certain genomic regions [34, 41, 70, 71].

The presence of a CTCF motif closely correlates with CTCF binding: over 90% of cohesin-binding sites overlapping with CTCF peaks contain CTCF motifs (Fig. 2c). In contrast, less than half of cohesin-binding sites harbor CTCF motifs in the absence of CTCF binding. Cohesin-binding sites without CTCF binding tend to be highly deviated from a CTCF motif, reflecting a CTCF-independent mechanism of recruitment (Fig. 2d). Interestingly, a small population of cohesin-CTCF overlapped sites also lack any CTCF motif, suggesting an alternative way by which cohesin and CTCF bind to these regions (Fig. 2c, d).

Nipbl reduction affects cohesin binding at CTCF-bound sites and repeat regions

In mESCs, it was proposed that Nipbl and CTCF recruit cohesin to different genomic regions, implying that cohesin binding to CTCF sites may be Nipbl-independent [40]. We noticed that when we ranked cohesin-binding sites based on the read number in wild type peaks, they matched closely with the ranking of cohesin-binding sites in mutant MEFs, indicating that the decrease of cohesin binding is roughly proportional to the strength of the original binding signals (Fig. 2e). This suggests that most cohesin-binding sites have similar sensitivity to Nipbl reduction. Importantly, CTCF-binding signals also correlate with the ranking of cohesin binding, indicating that CTCF-bound sites are in general better binding sites for cohesin (Fig. 2e). Because of this, they satisfy the peak definition despite the decrease of cohesin binding in mutant cells (Fig. 1f, g and Fig. 6b). This explains why CTCF-bound cohesin sites are apparently enriched in the sites that are common to both wild type and mutant cells (Fig. 2f).

Based on the above data, we further clarified the role of Nipbl in cohesin binding to CTCF sites. We compared the effect of Nipbl reduction on cohesin binding to representative sites, which have either CTCF binding or a CTCF motif or both (Fig. 3a). Decreased cohesin binding was observed at sites tested by manual ChIP-qPCR in Nipbl mutant MEFs, correlating with the decreased Nipbl binding (Fig. 3a). Consistent with the genome-wide ChIP-seq analysis (Fig. 1a), control histone H3 ChIP-qPCR revealed no significant differences at the corresponding regions, indicating that the decreased cohesin binding is not due to generally decreased ChIP efficiency in mutant MEFs compared to the wild type MEFs (Fig. 3a, bottom). Similar results were obtained using a small interfering RNA (siRNA) specific for Nipbl (Fig. 3c), which reduced Nipbl to a comparable level as in mutant cells (western blot in Fig. 3b and RT-qPCR results in Table 2). This demonstrates the specificity of the Nipbl antibody and confirms that the decreased cohesin binding seen in Nipbl mutant MEFs is the consequence of reduced Nipbl (Fig. 3a). Thus, Nipbl also functions in cohesin loading at CTCF sites.

Repeat sequences are often excluded from ChIP-seq analysis. However, cohesin binding is found at various repeat sequences, including pericentromeric and subtelomeric heterochromatin, and ribosomal DNA regions in the context of heterochromatin in mammalian cells [44, 45]. Thus, we also tested the effect of Nipbl reduction on cohesin binding to repeat sequences by manual ChIP-PCR (Fig. 3). Both Nipbl mutation (Fig. 3a, top) and Nipbl depletion by siRNA (Fig. 3b, c) resulted in decreased cohesin binding at the repeat regions, indicating that Nipbl is also important for cohesin binding to repeat sequences. In contrast, there were no significant differences in the histone H3 ChIP signals between these repeat regions in wild type and mutant MEFs (Fig. 3a, bottom). Taken together, the results indicate that Nipbl functions in cohesin loading even at CTCF sites and repeat regions, confirming the genome-wide decrease of cohesin binding caused by Nipbl haploinsufficiency.

Cohesin distribution patterns in the genome and enrichment in promoter regions

In order to gain insight into how the weakening of cohesin binding may affect gene expression in mutant cells, the distribution of cohesin-binding sites in the genomes of both wild type and mutant MEFs were examined. Approximately, 50% of all cohesin-binding sites are located in intergenic regions away from any known genes (Fig. 4a). However, there is a significant enrichment of cohesin binding in promoter regions, and to a lesser extent in the 3′ downstream regions, relative to the random genomic distribution generated by sampling from pre-immune ChIP-seq reads (Fig. 4b). Similar promoter and downstream enrichment has been observed in mouse and human cells [30, 31, 40, 42, 67] as well as in Drosophila [72]. Promoter enrichment is comparable in both wild type and Nipbl mutant MEFs, constituting ~ 10% of all the cohesin-binding sites (Fig. 4a). Thus, there is no significant redistribution or genomic region-biased loss of cohesin-binding sites in Nipbl mutant cells.

Cohesin-bound genes are sensitive to Nipbl haploinsufficiency

Based on the significant enrichment of cohesin binding in the promoter regions, we next examined the correlation between cohesin binding to the gene regions and the change of gene expression in mutant MEFs using a KS test. This is a nonparametric test for comparing peak binding sites with gene expression changes in the mutant MEFs (Fig. 5). Genes that displayed the greatest expression change in mutant MEFs compared to the wild type MEFs showed a strong correlation with cohesin binding to the gene region, indicating that direct binding to the target genes is the major mechanism by which cohesin mediates gene regulation in a Nipbl dosage-sensitive fashion (Fig. 5a, left). Random sampling of a comparable number of simulated peaks in the gene regions yielded no correlation (Fig. 5d, left). Interestingly, cohesin binding to the gene region correlates better with decreased gene expression than increased expression in mutant cells, indicating that gene activation, rather than repression, is the major mode of cohesin function at the gene regions (Fig. 5a, middle).

When analyzed separately, cohesin binding to the promoter regions (+ 2.5 kb to − 0.5 kb of transcription start sites (TSS) (Fig. 5a, right)) showed the highest correlation (p value = 3.3e−09) compared to the gene body and downstream (Fig. 5b). Thus, cohesin binding to the promoter regions is most critical for gene regulation. Similar to the entire gene region, cohesin binding correlates more significantly with a decrease in gene expression in mutant cells, which is particularly prominent at the promoter regions compared to gene bodies or downstream, indicating the significance of cohesin binding to the promoter regions in gene activation (Fig. 5c). Cohesin and CTCF binding closely overlapped at promoter regions in HeLa cells [31]. However, the overlap of CTCF binding with cohesin in MEFs is lower in the promoter regions (54%) than that in the intergenic regions (67%) [40]. Consistent with this, there is no significant correlation between CTCF binding in the promoter regions and gene expression changes in Nipbl mutant MEFs (p value = 0.28) by KS test (Fig. 5c, right). These results further indicate the cohesin-independent and Nipbl-insensitive function of CTCF in gene regulation. Taken together, the results suggest that cohesin binding to gene regions (in particular, to promoters) is significantly associated with gene activation that is sensitive to Nipbl haploinsufficiency.

Identification of cohesin target genes sensitive to Nipbl haploinsufficiency

Many of these Nipbl/cohesin-target genes contain cohesin-binding sites in more than one region (promoter, gene body and/or downstream), suggesting their collaborative effects (Fig. 6a). In particular, the promoter binding of cohesin is often accompanied by its binding to the gene body. However, binding pattern analysis revealed no significant correlation between a particular pattern and/or number of cohesin-binding sites and gene activation or repression (Fig. 6a). Rad21 ChIP-seq signal intensity profiles of several cohesin target genes (as defined above) reveal decreased cohesin binding in mutant cells at the binding sites originally observed in the wild type cells, supporting the notion that gene expression changes are the direct consequence of the reduced cohesin binding (Fig. 3a; Fig. 6b, top). There are other genes, however, that did not change expression significantly in mutant MEFs, but nevertheless also have reduced cohesin peaks nearby (Fig. 6b, bottom), suggesting that cohesin binding is not the sole determinant of the gene’s expression status and that its effect is context-dependent.

Nipbl- and cohesin-mediated activation of adipogenesis genes

One of the reported phenotypes of Nipbl+/− mice is their substantial reduction of body fat that mirrors what is observed in CdLS patients [19, 73]. It was found that Nipbl+/− MEFs exhibit dysregulated expression of several genes involved in adipocyte differentiation and reduced spontaneous adipocyte differentiation in vitro [19, 73]. We therefore examined the effect of Nipbl haploinsufficiency on these adipogenesis genes in detail. We found that many of them are bound by cohesin, in some cases at multiple sites, suggesting that cohesin plays a direct role in activation of these genes (Fig. 7). Although Il6 and Cebpδ were originally not included in the 115 genes due to low p values in the microarray analysis (Table 3 and Fig. 6a), significant expression changes were observed in mutant MEFs compared to the wild type MEFs by manual RT-qPCR. TNFα and PPARγ, also involved in adipogenesis, do not change their expression in mutant MEFs [19]. Importantly, a decrease of gene expression was observed not only in Nipbl+/− mutant cells but also by siRNA depletion of Nipbl, confirming that the effect is specifically caused by Nipbl reduction (Fig. 7a). Furthermore, depletion of cohesin itself decreased their expression even more significantly than Nipbl depletion. The results suggest that multiple genes involved in the adipogenesis pathway are direct cohesin targets that are sensitive to Nipbl haploinsufficiency.

Cohesin binding correlates significantly with H3K4me3 at the promoter

To investigate the genomic features associated with cohesin target genes, we examined the chromatin status of the target gene promoters. We found that cohesin peaks closely overlap with the peaks of H3K4me3, a hallmark of an active promoter, in a promoter-specific manner (Fig. 8a). In contrast, there are only minor peaks of H3K27me3 and even less H3K9me3 signal at cohesin-bound promoters. This is consistent with the results of the KS-test revealing the significant association of cohesin binding to the promoter regions with gene activation rather than repression (Fig. 5c). Interestingly, however, promoter binding of cohesin was found in genes with different expression levels in wild type MEFs, revealing no particular correlation with high gene expression (Fig. 8b). Cohesin target genes defined above (Table 3) also exhibit variable expression levels in wild type MEFs (Fig. 8b). Thus, their expression is altered in Nipbl mutant cells regardless of the original expression level in wild type cells, indicating that cohesin binding contributes to gene expression but does not determine the level of transcription per se.

When cohesin-bound genes were categorized in five different groups based on the gene expression status in wild type MEFs, significant H3K4me3 enrichment was observed even in the cohesin-bound promoters of genes with low expression, compared to cohesin-free promoters of genes with a similar expression level (Fig. 8c). Bivalent (H3K4me3 and H3K27me3) modifications are also enriched in the lowest gene expression category (Fig. 8c). Taken together, the results reveal that there is a close correlation between cohesin binding and H3K4me3 in the promoter regions regardless of the expression levels of the corresponding genes.

Reduced cohesin binding due to NIPBL reduction can lead to a loss of long-distance chromatin interaction

The above results revealed the critical association of cohesin binding to the promoter region and expression of the target genes. How does cohesin bound to the promoter affect gene expression? We recently showed that cohesin-mediated long-distance chromatin interaction between distal enhancer and promoter regions was reduced at the β-globin locus, resulting in reduced gene expression, in Nipbl mutant mice [35]. Thus, we tested the potential involvement of cohesin binding to the Cebpβ gene, one of the target adipogenesis genes described above, in such long-distance chromatin interaction(s) and whether it is affected by Nipbl reduction using chromosome conformation capture (3C) analysis (Fig. 9). We tested several flanking sites that are positive for cohesin and RNA polymerase II (pol II) binding as well as H3K4me1 and H3K4me3, the hallmarks for enhancers [74–76] (Fig. 9A). We observed that the Cebpβ promoter interacts with one such region (Fig. 9A, B, the site “c”). Although the site c is associated with only a weak Rad21 ChIP-seq signal, SMC1 and SMC3 ChIP-seq signals were found at the same region [67], confirming that this is an authentic cohesin-binding site (Fig. 9A). The results indicate a selectivity of chromatin interactions among neighboring cohesin-binding sites, revealing that not all proximal cohesin-binding sites interact with each other. Since the other two regions are also bound by CTCF, this may be due to the directionality of CTCF/cohesin binding [77, 78]. Importantly, the observed interaction is indeed reduced in both Nipbl mutant and Nipbl siRNA-treated MEFs (Fig. 9B). The 3C signals at the Cebpβ locus were normalized to the constant interaction observed at the Ercc3 locus [63, 64], which was not affected by Nipbl reduction. The results indicate that the decrease of long-distance chromatin interaction involving the promoters and distant DNA elements is one of the direct consequences of reduced cohesin binding, which may be one mechanism of gene expression alteration by Nipbl haploinsufficiency.

Nipbl functions in cohesin loading at both CTCF and non-CTCF sites

In mESCs, it was suggested that Nipbl is involved in cohesin binding to only a subset of cohesin-binding sites, which are largely distinct from CTCF-bound sites [40]. However, we found that Nipbl binds to, and its haploinsufficiency decreased cohesin binding to, CTCF sites in MEFs. A similar decrease of cohesin binding was observed at both CTCF insulators and non-CTCF sites in the β-globin locus in Nipbl+/− fetal mouse liver [35]. Furthermore, during differentiation in mouse erythroleukemia cells, both Nipbl and cohesin binding is concomitantly increased at these sites [35]. Therefore, while cohesin was suggested to slide from the Scc2 (Nipbl homolog)-dependent loading sites in yeast [79, 80], Nipbl is present and appears to directly affect cohesin loading at CTCF sites in mammalian cells. Nipbl, rather than cohesin, interacts with Mediator and HP1 and appears to recruit and load cohesin onto genomic regions enriched for Mediator and HP1 for gene activation and heterochromatin assembly, respectively [40, 45]. In contrast, cohesin, and not Nipbl, primarily interacts with CTCF [45, 81]. Thus, for cohesin binding to CTCF sites, we envision that cohesin initially recruits Nipbl that in turn stably loads cohesin onto CTCF sites.

A recent study indicated that almost all CTCF sites are bound by cohesin in primary mouse liver [43]. In MEFs, however, we found that ~ 42% of CTCF-bound sites appear to be cohesin-free. Furthermore, there is less overlap of cohesin and CTCF in the promoter regions compared to the intergenic regions, and little correlation between CTCF binding to the promoter and gene expression changes in Nipbl mutant cells was observed. Thus, in contrast to the cooperative function of cohesin and CTCF at distantly located insulator sites [36], cohesin and CTCF appear to have distinct functions at gene promoters. Distinct gene regulatory functions of CTCF and cohesin have also been reported in human cells [41]. Further study is needed to understand the recruitment specificity and functional relationship of cohesin and CTCF in gene regulation.

How does Nipbl haploinsufficiency affect cohesin target gene expression?

One mechanism of cohesin action in gene regulation is to mediate chromatin loop formation [35, 40]. Increased Nipbl and cohesin binding correlates with the induction of the enhancer-promoter interaction and robust gene activation at the β-globin locus [35]. Depletion of cohesin resulted in decreased enhancer-promoter interactions and downregulation of globin genes [35]. Similarly, Nipbl haploinsufficiency results in less cohesin binding and decreased promoter-enhancer interactions and β-globin gene expression [35]. In the current study, we also found that the cohesin-bound promoter of one of the target genes, Cebpβ, is involved in a long-distance chromatin interaction with a putative enhancer, which is decreased in Nipbl mutant cells, consistent with the decreased gene expression. Thus, Nipbl haploinsufficiency affects cohesin target gene expression by decreasing cohesin-mediated chromatin interactions.

It should be noted, however, that not all genes that we examined showed significant long-distance chromatin interactions involving cohesin-bound promoters. While this may be because we did not test the correct enhancer regions, it also suggests that cohesin may promote gene activation by a mechanism(s) other than by mediating long-distance promoter interaction. One possibility is gene looping. In Saccharomyces cerevisiae, the promoter and terminator regions of genes interact with each other, which was thought to facilitate transcription re-initiation [82]. Although cohesin is often found at the promoter and terminator regions of genes in MEFs, we failed to obtain any evidence for the involvement of these sites in gene looping with our limited analysis. Thus, how (or whether) cohesin at the promoter may regulate gene transcription in a loop formation-independent manner is currently unclear.

Cohesin binding to the gene body regions is found at many of the cohesin target genes. This may represent the cohesin binding at intragenic enhancer elements or may be related to Pol II pausing [29]. While cohesin was shown to facilitate Pol II elongation in Drosophila [83–85], cohesin together with CTCF in the intragenic region was found to cause Pol II pausing at the PUMA gene in human cells [86], suggesting that cohesin can have both positive and negative effects on transcriptional elongation in a context-dependent manner. Furthermore, not all the cohesin-bound genes changed expression in Nipbl+/− MEFs, echoing this notion that the effect of cohesin binding on gene expression is context-dependent. What determines the effects of cohesin binding at individual binding sites on gene expression requires further investigation.

The role of cohesin in the maintenance of gene expression

While there is now strong evidence for cohesin’s role in chromatin organization and gene activation, whether cohesin is involved in initiation or maintenance of gene activation is less clear. Enrichment of cohesin binding at the transcription start sites and termination sites was observed previously in mouse immune cells with no significant correlation to gene expression [30]. Our genome-wide analysis also revealed that cohesin binding to the gene regions has no obvious relationship to the level of gene expression in wild type MEFs. And yet, a decrease in cohesin binding is associated with a tendency to downregulate these genes, indicative of the positive role of cohesin on gene expression, consistent with the enriched presence of H3K4me3 in promoter regions. We speculate that cohesin may not be the primary determinant of gene activation, but rather cohesin binding may be important for maintaining gene expression status initially determined by sequence- and cell type-specific transcription factors. Similarly, enrichment of bivalent histone modifications in the promoters of cohesin-bound genes with very low expression suggests that cohesin also contributes to the maintenance of the poised state of these genes.

Nipbl haploinsufficiency vs. cohesin mutation

There are two different cohesin complexes in mammalian somatic cells that differ by one non-SMC subunit (i.e., SA1 (STAG1) or SA2 (STAG2)) [87, 88]. A recent report on SA1 knockout mice revealed some phenotypic similarity to what is seen in mice with Nipbl haploinsufficiency [67]. Interestingly, the SA1 gene is one of the cohesin target genes that is slightly upregulated in Nipbl mutant cells [19]. Thus, together with the compensatory increase of Nipbl expression from the intact allele, there appears to be a feedback mechanism that attempts to balance the expression of Nipbl and cohesin in response to Nipbl mutation. The fact that upregulation was observed with the SA1, but not SA2, gene may reflect the unique transcriptional role of SA1 [67]. Interestingly, however, only 10% of genes altered in Nipbl mutant MEFs are changed significantly in SA1 KO MEFs [67]. This discrepancy may, as observed in Drosophila [89], reflect the different effects of decreased binding versus complete knockout of a cohesin subunit on target gene expression. It could also be a result of the decreased binding of the second cohesin complex, cohesin-SA2.

Cohesin binding was relatively uniformly decreased genome-wide in Nipbl haploinsufficient cells with no significant redistribution of cohesin-binding sites. Point mutations of different subunits of cohesin cause CdLS and CdLS-like disorders with both overlapping and distinct phenotypes compared to CdLS cases caused by NIPBL mutations [9, 10, 13]. Non-overlapping effects of downregulation of different cohesin subunits have been reported in zebrafish [20, 26]. This may reflect an unequal role of each cohesin subunit in gene regulation, and it is possible that some of the cohesin target genes may be particularly sensitive to a specific cohesin subunit mutation. For example, similar to the TBP-associating factors (TAFs) in TFIID [90], cohesin subunits may provide different interaction surfaces for distinct transcription factors, which would dictate their differential recruitment and/or transcriptional activities. Furthermore, recent studies provide evidence for cohesin-independent roles of NIPBL in chromatin compaction and gene regulation [27, 28, 91]. Thus, disturbance of cohesin functions as well as impairment of cohesin-independent roles of NIPBL may collectively contribute to CdLS caused by NIPBL mutations.