Skip to main content
Download PDF

p63 cooperates with CTCF to modulate chromatin architecture in skin keratinocytes

Epigenetics & Chromatin volume 12, Article number: 31 (2019) Cite this article

Abstract

The transcription factor p63 regulates epidermal genes and the enhancer landscape in skin keratinocytes. Its molecular function in controlling the chromatin structure is, however, not yet completely understood. Here, we integrated multi-omics profiles, including the transcriptome, transcription factor DNA-binding and chromatin accessibility, in skin keratinocytes isolated from EEC syndrome patients carrying p63 mutations, to examine the role of p63 in shaping the chromatin architecture. We found decreased chromatin accessibility in p63- and CTCF-bound open chromatin regions that potentially contributed to gene deregulation in mutant keratinocytes. Cooperation of p63 and CTCF seemed to assist chromatin interactions between p63-bound enhancers and gene promoters in skin keratinocytes. Our study suggests an intriguing model where cell type-specific transcription factors such as p63 cooperate with the genome organizer CTCF in the three-dimensional chromatin space to regulate the transcription program important for the proper cell identity.

Introduction

Skin development and homeostasis requires tightly regulated epidermal keratinocyte proliferation, differentiation and apoptosis, and these processes are governed by cooperation of cis-regulatory elements, transcription factors (TFs), chromatin accessibility as well as higher-order chromatin organization [1,2,3,4,5]. The TF p63 is a key regulator of epidermal development. At the molecular level, p63 regulates a large number of genes, for example, p21 (CDKN1A) in cell cycle arrest [6], Fras1 in maintaining basement membrane integrity [7] and key genes such as keratins, filaggrin and loricrin required for epidermal morphogenesis and differentiation [1, 8]. Furthermore, p63 directly regulates chromatin factors: Satb1, Lsh and Brg1, to control the chromatin remodeling during epidermis development [9,10,11]. Recent studies showed that p63 exerts a crucial role in establishing the enhancer landscape [4, 5, 12, 13]. Through active enhancers, p63 cooperates with its co-regulating TFs to modulate transcriptional program required for epidermal homeostasis [4].

In human, heterozygous mutations of TP63 encoding p63 cause a spectrum of ectoderm-related disorders [14]. For example, ectrodactyly–ectodermal dysplasia–cleft lip/palate (EEC) syndrome is caused by point mutations located in the p63 DNA-binding domain and manifests ectodermal dysplasia with defects in the epidermis and epidermal-related appendages, limb malformation and cleft lip/palate. Five hotspot mutations affecting amino acids, R204, R227, R279, R280 and R304, cover approximately 90% of the all EEC syndrome cases [15]. Our previous study showed that mutant p63 resulted in a genome-wide redistribution of enhancers in keratinocytes established from EEC patients [13]. Consistently, the gene network analysis identified a significant co-expression gene module of ‘nucleosome assembly’, implying a less-organized chromatin structure in EEC syndrome keratinocytes [13]. How mutant p63 affects the chromatin structure is, however, not yet clear.

In this study, we characterized the chromatin accessibility using ATAC-seq in keratinocytes established from EEC patients carrying p63 mutations, in comparison with control keratinocytes. A clear difference in chromatin accessibility that correlated with the transcriptional dynamics was detected. Unexpectedly, strong enrichment of CTCF binding sites was observed at control-specific open chromatin regions. By combining published promoter capture Hi-C seq data, we found that CTCF and p63 were cooperatively involved in DNA loops to regulate epidermal genes. Our findings provide new insights into the coordinated regulatory role of CTCF and p63 in chromatin interactions in epidermal keratinocytes.

Results

Differential chromatin accessibility between control and p63 mutant keratinocytes

We performed assay for transposase-accessible chromatin with high-throughput sequencing (ATAC-seq) to characterize the accessible genome as well as the nucleosome position in both control and p63 mutant keratinocytes at the proliferation stage. Two replicas of ATAC-seq analyses showed high correlation (Fig. 1a), indicating high reproducibility. The principle component analysis (PCA) plot based on the chromatin accessibility displayed a clear separation between control and p63 mutant keratinocytes (Fig. 1b), which is highly consistent with results from gene expression profiles (Fig. 1b). The two p63 mutant lines showed high similarity in both genome accessibility and transcription level, as they are close on the PC1 axis that represents major variations. As the goal of this work is to examine the difference between control and EEC mutant lines, these two mutant lines were considered as one group, termed as p63 mutant keratinocytes, in the following analyses.

Fig. 1
figure 1

Differential chromatin accessibility in p63 mutant keratinocytes. a A heatmap of sample correlation matrix showing high similarities between duplicates and dissimilarities between control (CTR) and p63 mutant keratinocytes (R204W and R304W). b PCA plots of ATAC-seq and RNA-seq of control and p63 mutant keratinocytes at the proliferation stage. c Upper panel, volcano plots of ATAC-seq comparisons between control and p63 mutant keratinocytes. The x axis shows the log2 fold change of reads detected at the open chromatin regions (OCRs) and the y axis shows − log10 (p value). Lower panel, the transcriptional changes of differentially expressed genes associated with differential OCRs indicated in the volcano plots. d GO annotation of differentially expressed genes nearby Ctr-OCRs and Mt-OCRs. e Examples of Ctr-OCRs and their associated DE gene, CDC7, as well as Mt-OCRs and their associated DE gene, FAP. ATAC-seq tracks are shown in the RPKM scale

Full size image

Subsequently, we identified the differential open chromatin regions (OCRs) marked by ATAC-seq signals between control and p63 mutant keratinocytes. In total, there were 2492 open chromatin regions that showed higher signal in control keratinocytes, termed as control-specific open chromatin regions (Ctr-OCRs); in parallel, there were 3716 regions that showed higher signal in mutant keratinocytes, termed as mutant-specific open chromatin regions (Mt-OCRs) (Fig. 1c). As expected, we found that differential OCRs were positively correlated with gene expression, when they were assigned to the nearest differentially expressed (DE) genes (Fig. 1c). Genes associated with Ctr-OCRs were mainly involved in the regulation of ‘cell cycle’ (Fig. 1d), e.g., CDC7 (Fig. 1e), while genes nearby Mt-OCRs were mainly enriched in ‘actin filament-based process’ (Fig. 1d), e.g., FAP (Fig. 1e).

p63 and CTCF occupancy at Ctr-OCRs

To identify the potential TFs involved in these differential OCRs, we first performed a comparative motif analysis between control and p63 mutant keratinocytes using HOMER. Among the Ctr-OCRs, we found that p63 motif is most enriched (Fig. 2a), which is consistent with our previous finding that loss of p63 binding resulted in decreased enhancers in p63 mutant keratinocytes [13]. Unexpectedly, the motif of CTCF that is a well-known genome organizer [16] was identified as the second highly enriched motif in the Ctr-OCRs. Among the Mt-OCRs, AP-1 motif family was most enriched (Fig. 2b), consistent with findings on gained enhancers in p63 mutant keratinocytes [13].

Fig. 2
figure 2

Decreased chromatin accessibility and altered nucleosome organization at p63-bound OCRs and CTCF-bound OCRs in p63 mutant keratinocytes. a Top enriched motifs in Ctr-OCRs. 2492 Ctr-OCRs (targets) were used in this analysis with Mt-OCRs as background (bg). b Top enriched motifs in Mt-OCRs. 3716 Mt-OCRs (targets) were used in this analysis with Ctr-OCRs as background (bg). c Bandplots (top) and heatmaps (bottom) showing the quantification of ATAC-seq and p63 ChIP-seq signals at the 1223 p63-bound Ctr-OCRs. d Bandplots (top) and heatmaps (bottom) showing the quantification of ATAC-seq and CTCF ChIP-seq signals at the 376 CTCF-bound Ctr-OCRs. e V-plot of ATAC-seq reads showing the read length and distance from either p63 motif (top) or CTCF motif (bottom) in control and p63 mutant keratinocytes. The short fragments (< 150 bp) were indicative of nucleosome free DNA

Full size image

To validate our motif analysis result and explore the role of p63 in controlling Ctr-OCRs, we used published p63 ChIP-seq data from our previous work [13]. In control keratinocytes, a total of 1223 p63 binding sites (BSs) were covered by the 2492 Ctr-OCRs, showing a significant overlap (p < 0.001). When ATAC-seq and p63 ChIP-seq reads were mapped to the 1223 p63-bound Ctr-OCRs, we found a decrease in both chromatin accessibility and p63 binding signals at the p63-bound Ctr-OCRs in p63 mutant keratinocytes (Fig. 2c).

To confirm the role of CTCF in Ctr-OCRs, we performed CTCF ChIP-seq in both control and p63 mutant keratinocytes. In total, we found 12,116 CTCF BSs, among which only 363 were differential in CTCF binding between control and p63 mutant keratinocytes (p < 0.05), indicating that CTCF binding is rather stable. In addition, among all putative CTCF BSs, we found 376 overlapped with the 2492 Ctr-OCRs, which is statistically significant (p < 0.001). We then mapped ATAC-seq and CTCF ChIP-seq reads to the 376 CTCF-bound Ctr-OCRs. As expected, decrease in ATAC-seq signals was observed at the CTCF-bound Ctr-OCRs in p63 mutant keratinocytes. However, unexpectedly, there was no significant change in CTCF binding signals at the CTCF-bound Ctr-OCRs in p63 mutant keratinocytes, as compared to the control keratinocytes (Fig. 2d). Another unexpected finding was that there is few overlap (n = 435) out of 45,350 p63 BSs and 12,116 CTCF BSs (Additional file 1: Fig. S1A), although both motifs were identified in Ctr-OCRs. Consistently, there was no visible CTCF binding signal at all p63 BSs, and vice versa, no visible p63 binding at the CTCF BSs (Additional file 1: Fig. S1B). We also observed decreased p63 binding and largely unchanged CTCF binding signal in p63 mutant keratinocytes at the co-localized binding sites (Additional file 1: Fig. S1C).

Subsequently, we analyzed the nucleosome positioning at p63-bound and CTCF-bound Ctr-OCRs by mapping the distribution of ATAC-seq fragments centered on either p63 motif or CTCF motif (Fig. 2e). In control keratinocytes, we found an enrichment of short fragments (< 150 bp) surrounding both p63 motif and CTCF motif. The enrichment of such short reads clearly decreased in p63 mutant keratinocytes (Fig. 2e and Additional file 2: Fig. S2). We also observed an enrichment of fragments that had lengths slightly shorter than 200 bp which represent the stable flanking nucleosomes in control keratinocytes. However, in p63 mutant keratinocytes, the approximate 200 bp length fragments seemed to be less present (Additional file 2: Fig. S2), indicating an altered nucleosome organization.

Taken together, we observed decreased accessibility at p63-bound and CTCF-bound Ctr-OCRs in p63 mutant keratinocytes when compared to control keratinocytes. Notably, there was a clear decrease in p63 binding at the p63-bound Ctr-OCRs, whereas CTCF binding did not seem to be affected at the CTCF-bound Ctr-OCRs in p63 mutant keratinocytes.

Characterization of p63-bound OCRs and CTCF-bound OCRs

We reasoned that the chromatin states and the genomic locations of these Ctr-OCRs might be informative to understand the difference between p63-bound and CTCF-bound Ctr-OCRs. To this end, we firstly analyzed the chromatin states (CSs) of all CTCF and p63 BSs with the combination of six histone marks from Roadmap project [17] (Fig. 3a). We found that p63 bound to both promoters and enhancers marked by H3K27ac and H3K4me1, whereas CTCF was more enriched at promoters marked by H3K4me3 (Fig. 3a and Additional file 3: Fig. S3). We then performed analyses on the p63- and CTCF-bound Ctr-OCRs (Fig. 3b). When comparing to all Ctr-OCRs, p63-bound Ctr-OCRs were overrepresented at enhancers (CSs 8, 9 and 11), while the majority of CTCF-bound Ctr-OCRs showed overrepresentation at promoter-proximal regions (CSs 2 and 4, Fig. 3b), consistent with all p63 and CTCF BSs.

Fig. 3
figure 3

Chromatin state and locations of p63- and CTCF-bound OCRs. a Fold enrichment/depletion of chromatin states (CSs) in p63 BSs and CTCF BSs compared to the randomized sequences. The 18 CSs were defined using ENCODE data [17]. b Relative enrichment (> 1)/depletion (< 1) of CSs in p63-bound Ctr-OCRs and CTCF-bound Ctr-OCRs compared to all Ctr-OCRs. The circle size represents the percentage of either p63-bound Ctr-OCRs (left) or CTCF-bound Ctr-OCRs (right) in each CS. Relative score is calculated by comparing odds ratios between two overlap tests, in which odds ratio is estimated by BEDtools fisher option. c Percentage of different loop types defined in the PCHiC data. d Size distribution of promoter–enhancer loops in the PCHiC data. e Estimated enrichment of CTCF and p63 at promoters, enhancers, or TAD Boundary compared to the reshuffled sequences (CTR). Odds ratio indicates the enrichment (> 1)/depletion (< 1) of two independent peak lists calculated by BEDtools fisher option

Full size image

Next, to analyze the localization of p63-bound and CTCF-bound Ctr-OCRs in regulatory chromatin interactions, we utilized the promoter capture Hi-C (PCHiC) data and TAD (topologically associating domain) regions from published Hi-C data in keratinocytes [18]. In the PCHiC, promoters were used as baits to capture the interacting genomic regions. In total, this PCHiC dataset captured 119,648 interactions, including 103,497 (86.5%) promoter–enhancer loops and 16,151 (13.5%) promoter–promoter loops (Fig. 3c). As p63 plays a more prominent role at enhancers (Fig. 3a) [4], we focused on promoter–enhancer loops. These loops had a median loop size of 250.9 kb (Fig. 3d). Among the previously defined 1223 p63-bound Ctr-OCRs (Fig. 2c), 388 of them were located at the anchors of 1674 loops, at either the promoter or the enhancer, and on average one p63-bound Ctr-OCR was connected by ~ 4.3 loops. Similarly, among the 376 CTCF-bound Ctr-OCRs (Fig. 2d), 191 of them were located at the anchors of 1093 loops; on average, one CTCF-bound Ctr-OCR was connected by ~ 5.7 loops. Therefore, CTCF-bound Ctr-OCRs at anchors were connected by more loops.

Furthermore, we found that CTCF-bound Ctr-OCRs located at loop anchors were more enriched for promoters (fourfold) (Fig. 3e), consistent with all CTCF BSs and all CTCF Ctr-OCRs. In contrast, p63-bound Ctr-OCRs associated with loops were enriched in both promoters and enhancers to a similar extent (Fig. 3e), which is different from all p63 BSs and all p63-bound Ctr-OCRs that were only enriched at enhancers (Fig. 3b). We also examined whether p63- and CTCF-bound Ctr-OCRs were localized at the TAD boundaries. As expected, CTCF-bound Ctr-OCRs showed ~ 2.5-fold enrichment for TAD boundaries, whereas p63-bound Ctr-OCRs were depleted of TAD boundaries (Fig. 3e). Randomization by reshuffling genomic regions with same sizes did not show any enrichment in TAD boundaries, promoters or enhancers. As CTCF-bound Ctr-OCRs were enriched at the promoters and p63-bound Ctr-OCRs were enriched at both promoters and enhancers, it is plausible that p63-bound enhancers interact with CTCF-bound promoters to regulate genes through DNA looping and a subset of these loops may be affected by p63 mutations.

p63 and CTCF mediate a subset of loops to regulate transcription

To assess how p63 and CTCF cooperates via DNA looping to regulate gene expression, we performed further analyses using previously published gene expression data [13] to examine the relation between gene expression and DNA looping. In control keratinocytes, we defined highly expressed genes by quartiles to four groups using their transcriptional levels (Q1-4, FPKM ≥ 1), and another group with lowly expressed genes (Q0, FPKM < 1). We detected an increasing number of associated loops with elevated gene expression levels (Fig. 4a). This suggests that highly expressed genes may be more indispensable of regulation through loops, which is line with previous study showing multi-loop activation hubs at key regulatory genes [19]. Moreover, the differentially expressed genes between control and p63 mutant keratinocytes (Fig. 1b) were associated with significantly more loops (p < 0.001) than all expressed genes (Fig. 4b). We also performed the same analysis among loops with either p63-bound or CTCF-bound Ctr-OCRs. Similar results were found, although in these comparisons the difference was not significant, probably due to the low number of associated loops (Fig. 4c, d).

Fig. 4
figure 4

p63- and CTCF-bound Ctr-OCRs are involved in chromatin loops regulating epidermal genes. a Highly expressed genes are associated with more DNA loops. Genes were divided in quarters according to their expression levels in keratinocytes. The number of connected DNA loops associated with each group of genes were plotted in the bottom bar chart. b DE genes were associated with more connected DNA loops when compared to all expressed genes. ***, p < 0.001, Wilcoxon signed-rank test. c DE genes were associated with more connected DNA loops with p63-bound Ctr-OCRs when compared to all expressed genes. d DE genes were associated with more connected DNA loops with CTCF-bound Ctr-OCRs when compared to all expressed genes. e Overlapping of genes mapped by CTCF-bound Ctr-OCRs (orange) and p63-bound Ctr-OCRs (blue). f GO annotation of the 96 genes potentially co-regulated by p63 and CTCF. g p63-bound and CTCF-bound Ctr-OCRs were involved in the same or different loops connecting the promoters of the 96 genes. 60 genes were regulated by the same loops with both p63-bound and CTCF-bound Ctr-OCRs: (1) in 12.6% cases, p63-bound Ctr-OCR located at one anchor while CTCF-bound Ctr-OCR located at the other anchor of the loop; in the majority cases, p63-bound and CTCF-bound Ctr-OCRs were both present at the (2) promoter (70.6%) or (3) enhancer (16.8%) of the DNA loops. h The left 36 genes were co-regulated by connecting DNA loops with p63-bound and CTCF-bound Ctr-OCRs at enhancers of different loops. i PIGV is shown as an example of G (1), which is connected by the DNA loop pair1402, bridging the PIGV promoter (with a CTCF-bound Ctr-OCR) and a distal enhancer (with a p63-bound Ctr-OCR). j KRT5 is shown as an example of G (3), which is connected by the DNA loop pair25496, bridging the KRT5 promoter and a distal enhancer (with both CTCF-bound and p63-bound Ctr-OCRs). The blue arrows pointed to the p63 BSs while the green arrows pointed to the CTCF BSs

Full size image

We then mapped 388 p63-bound and 191 CTCF-bound Ctr-OCRs which are located at loop anchors to their target genes by identifying their interacting gene promoters located at anchors of the other end of the loops. In some cases, Ctr-OCRs were located in promoters, and the associated genes were considered as target genes. In total, 388 p63-bound Ctr-OCRs interacted with 304 gene promoters, while 191 CTCF-bound Ctr-OCRs interacted with 500 gene promoters through looping interactions, indicating that on average a CTCF-bound Ctr-OCR is involved in the regulation of more genes than a p63-bound Ctr-OCR. There was a significant (p < 0.001) overlap (n = 96) of genes that were potentially regulated by p63- and CTCF-bound Ctr-OCRs (Fig. 4e). These 96 genes represent genes connected by the promoter–enhancer loops that were anchored by both p63-bound and CTCF-bound Ctr-OCRs and that may be affected in p63 mutant keratinocytes. As the three-dimensional chromatin landscape is relatively stable, we reasoned that the difference of the chromatin landscape in control and p63 mutant keratinocytes may influence gene expression during keratinocyte differentiation [18]. Therefore, we examined the expression of these 96 genes using previously reported differential gene expression in control and p63 mutant keratinocytes during differentiation [13]. Indeed, we found that 39 of these 96 genes showed deregulation (p < 0.001). Notably, many of the 96 genes were involved in ‘epidermal cell differentiation,’ e.g., KRT5, KRT75 and B9D1 (Fig. 4f).

Next, we dissected how p63-bound and CTCF-bound Ctr-OCRs at loop anchors co-regulated these 96 genes: whether p63-bound and CTCF-bound Ctr-OCRs were involved in the same or different promoter–enhancer loops. We found that 60 genes were regulated by promoter–enhancer loops where p63-bound and CTCF-bound Ctr-OCRs were both located at the anchors in the same loops, although they may co-localize in different manner (Fig. 4g). Among these 60 genes, 26 showed differential expression between control and p63 mutant keratinocytes during differentiation, with PIGV and KRT5 shown as examples (Fig. 4i, j). Among the loops associated with the 60 genes, the case where p63-bound Ctr-OCR was located at one anchor, while CTCF-bound Ctr-OCR was located at the other anchor of the loop counted for 12.6% (Fig. 4g, (1)). In the majority of the cases, both p63-bound and CTCF-bound Ctr-OCRs were present either at the promoter (70.6%) (Fig. 4g, (2)) or at the enhancer (16.8%) (Fig. 4g, (3)). Interestingly, for the rest 36 genes gene promoters were connected with both p63-bound and CTCF-bound Ctr-OCRs from different loops (Fig. 4h), suggesting that enhancers occupied by p63-bound and CTCF-bound Ctr-OCRs can co-regulate the same gene through different loops. Among these 36 genes, 13 were deregulated in p63 mutant keratinocytes during differentiation.

Discussion

Over the years, molecular understanding of p63 function in keratinocytes has led to the identification of numerous target genes including chromatin factors [20,21,22] as well as of a solid role of p63 in orchestrating the enhancer landscape [4, 5, 23]. However, whether p63 is involved in the direct regulation of chromatin architecture is not known. In this study, we characterized the chromatin accessibility with ATAC-seq in both control and EEC syndrome patient keratinocytes carrying p63 mutations (R204W and R304W). Interestingly, we found decreased chromatin accessibility in p63- and CTCF-bound open chromatin regions in mutant p63 keratinocytes. These less accessible chromatin regions interact with epidermal gene promoters and potentially regulate gene expression. Collectively, we proposed an unreported role of p63 cooperating with CTCF in chromatin interactions.

It is envisaged that p63 is a key factor in skin keratinocytes [4, 5, 24]. Ample studies have shown that the proper regulation of gene expression by p63 through the precise control of enhancers is essential for maintaining the epidermal cell identity [4, 13, 25, 26]. It is equally important that p63 requires additional co-regulating TFs to regulate transcription. Several studies reported p63 co-regulating TFs, e.g., AP1, AP2, STAT5 and RFX5 [20, 22, 23, 27], which are essential for the epidermal gene expression in a temporal and spatial manner [4, 28]. Besides, p63 interacts with the chromatin remodeling factor BAF1 to maintain the open chromatin regions [5] and with Dnmt3a to locate enhancers [12]. In this study, we found that CTCF, a well-known TF maintaining the chromatin architecture, is also a p63 co-regulating TF modulating a subset of chromatin loops. However, this cooperation is probably not through direction protein–protein interactions but mainly through bridging enhancers and promoters, given the few overlap observed between p63 and CTCF BSs. It is known that the majority of p63-bound regions are located in intergenic or intragenic regions (enhancers) rather than promoters [22], and therefore, how p63-bound enhancers interact with gene promoters is an interesting question. Our current study provides the first line of evidence that the cooperation with CTCF assists the looping of p63-bound enhancers to gene promoters.

The observation of the enriched CTCF motif in Ctr-OCRs detected in our ATAC-seq analyses was somewhat unexpected. Our previous enhancer-centered analysis that focused on H3K27ac enriched genomic regions did not capture CTCF motifs [4, 13], probably due to the fact that CTCF is more enriched in open chromatin regions depleted of H3K27ac but enriched for H3K4me3 signal that are often promoters (Fig. 3a and Additional file 3: Fig. S3). Interestingly, detecting the CTCF motif in Ctr-OCRs in keratinocytes is in line with the previous finding that the motif of YY1 was enriched in p63-bound sequences [27]. YY1 has been shown to co-localize with CTCF and stabilize chromatin looping [29, 30].

CTCF maintains the genome organization by occupying the boundaries of megabase-scale TADs and functions as an insulator to block enhancer–promoter interactions between different TADs [16]. Our data showed that, in addition to demarcating TADs, CTCF mediates promoter–enhancer loops, often located in promoter-proximal regions (Fig. 3e), to facilitate the promoter–enhancer interactions within one TAD. This is in line with the concept that a subpopulation of CTCF associates with the RNA polymerase II (Pol II) protein complex to activate transcription [31]. It is likely that CTCF helps to bridge the p63-bound enhancers to transcription start site-proximal regulatory elements and to initiate transcription by interacting with Pol II, thus supporting a role of CTCF in facilitating contacts between transcription regulatory sequences. This model has been demonstrated by the previous work on the beta-globin locus [32]. To get a clearer picture, technologies that are similar to ChIA-PET [33], such as HiChIP [34], would be useful to understand long-range contacts associated with CTCF or p63 in skin keratinocytes.

The three-dimensional chromatin landscape is relatively stable once established in a specific cell type, probably during cell commitment, and cell type-specific looping structures mainly control the accessibility of enhancers to their specific targets [35]. Recent studies showed that cell type-specific TFs are involved in regulation of DNA looping in macrophage development, neural development as well as during cell reprogramming, which shed light on the role of cell type-specific TFs in regulating chromatin interactions [19, 36, 37]. Here, we proposed that a number of loci nearby epidermal genes were organized into a ‘regulatory chromatin hub’ within the chromatin interactions mediated by CTCF in epidermal keratinocytes (Fig. 5a). Such hubs contain multiple connecting DNA loops that require not only CTCF binding that is rather static but also binding of cell type-specific TFs (e.g., p63) for the transcriptional activity. In this model, cell type-specific TFs may be essential to make the DNA loops active in transcription. This hypothesis is consistent with our observation that unchanged CTCF binding signal but decreased accessibility in the CTCF-bound Ctr-OCRs in p63 mutant keratinocytes (Fig. 2d).

Fig. 5
figure 5

Models of how p63 and CTCF co-mediate chromatin loops regulating epidermal genes. a In control keratinocytes, p63 and CTCF co-mediate a subset of DNA loops nearby epidermal genes and regulate their gene expression. b In p63 mutant keratinocytes, due to the loss of p63 binding, p63 and CTCF co-mediated DNA loops are not functional anymore, and therefore epidermal genes are deregulated. There are two possibilities in deregulated function of looping: (1) the looping structure is stable as CTCF binding is unchanged; (2) the looping structure is partially changed because of loss of p63 binding and loss of the p63-dependent loop

Full size image

Moreover, our observations shed light on the chromatin-based mechanisms underlying EEC syndrome caused p63 mutations. These EEC mutations, such as R204W and R304W, were shown to disrupt p63 DNA-binding, resulting in impaired transactivation activity and loss of epidermal cell identity [13, 38, 39]. Our data suggest that deregulated function of DNA loops mediated by p63 and CTCF represents an additional layer to the disease mechanism (Fig. 5b). Given that CTCF binding did not change, it is likely that the looping interactions are static, but the looping is no longer active due to the absence of p63 binding (Fig. 5b (1)). Therefore, the epidermal genes were down-regulated and p63 mutant keratinocytes showed loss of epidermal cell identity [13]. On the other hand, it is also possible that the loop ablation as a result of lost p63 binding leads to the decreased accessibility at CTCF-bound promoters and deregulated gene expression in p63 mutant keratinocytes (Fig. 5b (2)). The latter scenario suggests that p63 is the key factor determining the regulatory activity of such promoter–enhancer interactions. Nevertheless, how p63 mutations affect the cooperation between p63 and CTCF and gene expression is still unclear. Among the 96 genes with loops connected by p63- and CTCF-bound Ctr-OCRs, only 39 of them were deregulated in p63 mutant keratinocytes. It is possible that ‘shadow enhancers’ that are composed of clusters of enhancer and contribute to robust gene transcription [40, 41], are involved in regulation of the 57 genes that are not deregulated in mutant keratinocytes. Further experiments are required to test these hypotheses.

Taken together, our results highlight the cooperation of p63 and CTCF in mediating DNA loops as well as gene regulation. These observations provide new insights into the general regulatory role of chromatin looping which depends on the coordination of CTCF and cell type-specific TFs such as p63 in skin keratinocytes.

Materials and methods

Ethics statement

All procedures for establishing and maintaining human primary keratinocytes were approved by the ethical committee of the Radboud university medical center (“Commissie Mensgebonden Onderzoek Arnhem-Nijmegen”). Informed consent was obtained from all donors of a skin biopsy.

Human primary keratinocyte culture

Primary keratinocytes were established previously from skin biopsies of three EEC syndrome patients carrying heterozygous mutations in the p63 DNA-binding domain, R204W [42] and R304W [39], as well as of non-EEC volunteers (Dombi23, referred to as control) [43]. R204W and R304W are treated as a p63 mutant group for further analyses to minimize the effect from individual difference. As previously described [44], primary keratinocytes were cultured in Keratinocyte Basal Medium (KBM, Lonza #CC-4131) supplemented with 100 U/mL penicillin/streptomycin (Gibco Life Technology #15140122), 0.1 mM ethanolamine (Sigma-Aldrich #141-43-5), 0.1 mM O-phosphoethanolamine (Sigma-Aldrich #1071-23-4), 0.4% (vol/vol) bovine pituitary extract, 0.5 μg/mL hydrocortisone, 5 μg/mL insulin and 10 ng/mL epidermal growth factor (Lonza #CC-4131). Medium was refreshed every other day. When cells were more than 90% confluent, cells were collected at the proliferation stage. No mycoplasma contamination is found during cell culture.

ATAC-seq

ATAC libraries of the control and p63 mutant keratinocytes were prepared by a documented protocol [45]. In brief, keratinocytes were treated with Accutase® solution on plates and well re-suspended into single cells. Pellet cells for 5 min at 500 g. After twice wash with ice-old PBS, pellet the cells again. Re-suspend cells and take out 100,000 cells for lysis with ice-old freshly made lysis buffer (containing 10 mM Tris–HCl pH 7.5, 10 mM NaCl, 3 mM MgCl2 and 0.1% IGEPAL CA-630 detergent). Then perform tagmentation using 2 μL of Tn5 transposase and 12.5 uL 2 × TD buffer (Illumina #FC-121-1031) at 37 °C for 1 h with 650 rpm shaking. The resulted DNA fragments underwent two sequential seven-cycle PCR amplifications, and in between the libraries were selected for < 500 bp fragments using SPRI beads. The final PCR products were purified with QIAquick PCR Purification Kit (QIAGEN #28106) and quantified with the KAPA Library Quantification Kit (Kapa Biosystems #KK4844), and then sequenced in a paired-ended manner using the NextSeq 500 (Illumina) according to standard Illumina protocols.

CTCF ChIP-seq

Chromatin for ChIP was prepared as previously described [20]. ChIP assays were performed following a standard protocol [46] with minor modifications. On average, 0.5 M keratinocytes were used in each ChIP. 4x ChIP reactions are pooled to prepare one ChIP-seq sample. Antibodies against CTCF (Millipore #07-729, 5uL) were used in each ChIP assay. Resulted DNA fragments from four independent ChIP assays were purified with QIAquick PCR Purification Kit (QIAGEN #28106). Afterward, 5 ng DNA fragments were pooled and proceeded on with library construction using KAPA Hyper Prep Kit (Kapa Biosystems #KK8504) according to the standard protocol. The prepared libraries were then sequenced using the NextSeq 500 (Illumina) according to standard Illumina protocols.

RNA-seq data analysis

Paired-end RNA-seq reads were mapped to the human genome hg19 using STAR aligner [47] in two-pass mode with default parameters and enumerate stranded gene-level read counts at the same time. The generated count matrix was used as input for DESeq 2 package [48] to distinguish differential expressed genes between control and mutant keratinocytes. These genes greater than 1.5-fold changed at adjusted p value < 0.05 were considered significantly deregulated.

ChIP-seq and ATAC-seq data analysis

Sequenced reads were aligned against the UCSC hg19 human reference genome with Burrows–Wheeler Aligner (BWA) program with default parameters [49]. The potential PCR and optical duplicates were removed using Picard MarkDuplicates option. The filtered BAM files were inputted to MACS2 [50] for peak calling. The p63 ChIP-seq data were published and re-analyzed in this study [13]. The p63, CTCF and H3K27ac were called using the narrow setting (default) with a q value of 0.01. For ATAC-seq, the open chromatin regions were predicted using default parameters except for using -f BAMPE option. Peaks overlapping with the consensus excludable ENCODE blacklist were dropped to avoid confounding by repetitive regions. All alignment files were extended to 200-bp and scaled to RPKM-normalized read coverage files using deepTools [51] for visualization. To compare binding profiles between different samples unbiasedly, we applied library size factors estimated from DESeq2 [48] on RPKM values.

Differentially accessible regions were detected using DESeq 2 package [48] with fold change less than 2.0 and p value below 0.05. Differential motif analysis in the differential DHSs was employed by the findMotifs function in HOMER tool (http://homer.salk.edu/homer/motif/) with other default parameters, which can normalize the background sequences to remove GC-bias. The BEDtools suite (https://bedtools.readthedocs.io/en/latest/content/bedtools-suite.html) was used to test overlap and enrichment between different intervals.

Dimensionality reduction and functional annotation

For visualization, top 1000 variable genes were first selected based on interquartile range (IQR) of normalized expression values, and further used to reduce dimensionality (principal component analysis) of the dataset by pca function in R. Functional enrichment was evaluated by Metascape online tool [52], to gain insight into the biological functions for deregulated genes. Only these functional terms with Benjamini-adjusted p value < 0.05 were considered significantly overrepresented.

Capture Hi-C data processing

Raw promoter capture Hi-C sequencing reads from GSE84662 were processed using the HiCUP pipeline [53], like quality control, alignment to hg19 and reads filtering. Technical replicates were merged and de-duplicated, and then, pooled biological replicates were then jointly used for interaction identification with CHiCAGO [54] and the associated chicagoTools suite. All significant interactions were defined based on a strict interaction threshold (CHiCAGO score ≥ 5).

Briefly, CHiCAGO predicts interactions based on a convolution background model reflecting both ‘Brownian’ (real, but expected interactions) and ‘technical’ (assay and sequencing artifacts) components. The putative p values are corrected using a weighted false discovery control procedure that specifically accommodates the fact that increasingly larger numbers of tests are performed at regions where progressively smaller numbers of interactions are expected. The weights were learned based on the decrease in the reproducibility of interaction calls between the individual replicates of macrophage samples with distance. Interaction scores were then computed for each fragment pair as log-transformed, soft-thresholded, weighted p values.

References

  1. Truong AB, Kretz M, Ridky TW, Kimmel R, Khavari PA. p63 regulates proliferation and differentiation of developmentally mature keratinocytes. Genes Dev. 2006;20(22):3185–97.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  2. Lippens S, Denecker G, Ovaere P, Vandenabeele P, Declercq W. Death penalty for keratinocytes: apoptosis versus cornification. Cell Death Differ. 2005;12(S2):1497.

    Article  CAS  PubMed  Google Scholar 

  3. Botchkarev VA, Flores ER. p53/p63/p73 in the epidermis in health and disease. Cold Spring Harbor Perspect Med. 2014;4(8):a015248.

    Article  Google Scholar 

  4. Kouwenhoven EN, Oti M, Niehues H, van Heeringen SJ, Schalkwijk J, Stunnenberg HG, van Bokhoven H, Zhou H. Transcription factor p63 bookmarks and regulates dynamic enhancers during epidermal differentiation. EMBO Rep. 2015;16(7):863–78.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Bao X, Rubin AJ, Qu K, Zhang J, Giresi PG, Chang HY, Khavari PA. A novel ATAC-seq approach reveals lineage-specific reinforcement of the open chromatin landscape via cooperation between BAF and p63. Genome Biol. 2015;16(1):284.

    Article  PubMed  PubMed Central  Google Scholar 

  6. LeBoeuf M, Terrell A, Trivedi S, Sinha S, Epstein JA, Olson EN, Morrisey EE, Millar SE. Hdac1 and Hdac2 act redundantly to control p63 and p53 functions in epidermal progenitor cells. Dev Cell. 2010;19(6):807–18.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. Koster MI, Dai D, Marinari B, Sano Y, Costanzo A, Karin M, Roop DR. p63 induces key target genes required for epidermal morphogenesis. Proc Natl Acad Sci. 2007;104(9):3255–60.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. Romano R-A, Ortt K, Birkaya B, Smalley K, Sinha S. An active role of the ΔN isoform of p63 in regulating basal keratin genes K5 and K14 and directing epidermal cell fate. PLoS ONE. 2009;4(5):e5623.

    Article  PubMed  PubMed Central  Google Scholar 

  9. Fessing MY, Mardaryev AN, Gdula MR, Sharov AA, Sharova TY, Rapisarda V, Gordon KB, Smorodchenko AD, Poterlowicz K, Ferone G, et al. p63 regulates Satb1 to control tissue-specific chromatin remodeling during development of the epidermis. J Cell Biol. 2011;194(6):825–39.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Mardaryev AN, Gdula MR, Yarker JL, Emelianov VU, Poterlowicz K, Sharov AA, Sharova TY, Scarpa JA, Joffe B, Solovei I, et al. p63 and Brg1 control developmentally regulated higher-order chromatin remodelling at the epidermal differentiation complex locus in epidermal progenitor cells. Development. 2014;141(1):101–11.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Keyes WM, Pecoraro M, Aranda V, Vernersson-Lindahl E, Li W, Vogel H, Guo X, Garcia EL, Michurina TV, Enikolopov G. ΔNp63α is an oncogene that targets chromatin remodeler Lsh to drive skin stem cell proliferation and tumorigenesis. Cell Stem Cell. 2011;8(2):164–76.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Rinaldi L, Datta D, Serrat J, Morey L, Solanas G, Avgustinova A, Blanco E, Pons JI, Matallanas D, Von Kriegsheim A, et al. Dnmt3a and Dnmt3b associate with enhancers to regulate human epidermal stem cell homeostasis. Cell Stem Cell. 2016;19(4):491–501.

    Article  CAS  PubMed  Google Scholar 

  13. Qu J, Tanis S, Smits J, Kouwenhoven EN, Oti M, Logie C, Stunnenberg H, Mulder K, Zhou H. Mutant p63 affects epidermal cell identity through rewiring the enhancer landscape. Cell Rep. 2018;25:3490–503.

    Article  CAS  PubMed  Google Scholar 

  14. Rinne T, Brunner HG, van Bokhoven H. p63-associated disorders. Cell Cycle. 2007;6(3):262–8.

    Article  CAS  PubMed  Google Scholar 

  15. Rinne T, Hamel B, Bokhoven HV, Brunner HG. Pattern of p63 mutations and their phenotypes—update. Am J Med Genet A. 2006;140(13):1396–406.

    Article  PubMed  Google Scholar 

  16. Ong C-T, Corces VG. CTCF: an architectural protein bridging genome topology and function. Nat Rev Genet. 2014;15(4):234.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Kundaje A, Meuleman W, Ernst J, Bilenky M, Yen A, Heravi-Moussavi A, Kheradpour P, Zhang Z, Wang J, Ziller MJ. Integrative analysis of 111 reference human epigenomes. Nature. 2015;518(7539):317.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Rubin AJ, Barajas BC, Furlan-Magaril M, Lopez-Pajares V, Mumbach MR, Howard I, Kim DS, Boxer LD, Cairns J, Spivakov M. Lineage-specific dynamic and pre-established enhancer–promoter contacts cooperate in terminal differentiation. Nat Genet. 2017;49(10):1522.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Phanstiel DH, Van Bortle K, Spacek D, Hess GT, Shamim MS, Machol I, Love MI, Aiden EL, Bassik MC, Snyder MP. Static and dynamic DNA loops form AP-1-bound activation hubs during macrophage development. Mol Cell. 2017;67(6):1037–48.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Kouwenhoven EN, van Heeringen SJ, Tena JJ, Oti M, Dutilh BE, Alonso ME, de la Calle-Mustienes E, Smeenk L, Rinne T, Parsaulian L, et al. Genome-wide profiling of p63 DNA-binding sites identifies an element that regulates gene expression during limb development in the 7q21 SHFM1 locus. PLoS Genet. 2010;6(8):e1001065.

    Article  PubMed  PubMed Central  Google Scholar 

  21. Vigano MA, Lamartine J, Testoni B, Merico D, Alotto D, Castagnoli C, Robert A, Candi E, Melino G, Gidrol X, et al. New p63 targets in keratinocytes identified by a genome-wide approach. EMBO J. 2006;25(21):5105–16.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. McDade SS, Henry AE, Pivato GP, Kozarewa I, Mitsopoulos C, Fenwick K, Assiotis I, Hakas J, Zvelebil M, Orr N, et al. Genome-wide analysis of p63 binding sites identifies AP-2 factors as co-regulators of epidermal differentiation. Nucleic Acids Res. 2012;40(15):7190–206.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Sethi I, Sinha S, Buck MJ. Role of chromatin and transcriptional co-regulators in mediating p63-genome interactions in keratinocytes. BMC Genom. 2014;15:1042.

    Article  Google Scholar 

  24. Cavazza A, Miccio A, Romano O, Petiti L, Malagoli Tagliazucchi G, Peano C, Severgnini M, Rizzi E, De Bellis G, Bicciato S, et al. Dynamic transcriptional and epigenetic regulation of human epidermal keratinocyte differentiation. Stem Cell Reports. 2016;6(4):618–32.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Fan X, Wang D, Burgmaier JE, Teng Y, Romano RA, Sinha S, Yi R. Single cell and open chromatin analysis reveals molecular origin of epidermal cells of the skin. Dev Cell. 2018;47:21–37.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Soares E, Xu Q, Li Q, Qu J, Zheng Y, Raeven HHM, Brandao K, van den Akker WMR, Tang F, Zhou H. Single-cell RNA-seq identifies a reversible epithelial-mesenchymal transition in abnormally specified epithelia of p63 EEC syndrome. bioRxiv 2019;437632.

  27. Yang A, Zhu Z, Kapranov P, McKeon F, Church GM, Gingeras TR, Struhl K. Relationships between p63 binding, DNA sequence, transcription activity, and biological function in human cells. Mol Cell. 2006;24(4):593–602.

    Article  CAS  PubMed  Google Scholar 

  28. Kouwenhoven EN, van Bokhoven H, Zhou H. Gene regulatory mechanisms orchestrated by p63 in epithelial development and related disorders. Biochim Biophys Acta. 2015;1849(6):590–600.

    Article  CAS  PubMed  Google Scholar 

  29. Schwalie PC, Ward MC, Cain CE, Faure AJ, Gilad Y, Odom DT, Flicek P. Co-binding by YY1 identifies the transcriptionally active, highly conserved set of CTCF-bound regions in primate genomes. Genome Biol. 2013;14(12):R148.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Beagan JA, Duong MT, Titus KR, Zhou L, Cao Z, Ma J, Lachanski CV, Gillis DR, Phillips-Cremins JE. YY1 and CTCF orchestrate a 3D chromatin looping switch during early neural lineage commitment. Genome Res. 2017;27:1139–52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Chernukhin I, Shamsuddin S, Kang SY, Bergström R, Kwon Y-W, Yu W, Whitehead J, Mukhopadhyay R, Docquier F, Farrar D. CTCF interacts with and recruits the largest subunit of RNA polymerase II to CTCF target sites genome-wide. Mol Cell Biol. 2007;27(5):1631–48.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Wendt KS, Grosveld FG. Transcription in the context of the 3D nucleus. Curr Opin Genet Dev. 2014;25:62–7.

    Article  CAS  PubMed  Google Scholar 

  33. Fullwood MJ, Wei C-L, Liu ET, Ruan Y. Next-generation DNA sequencing of paired-end tags (PET) for transcriptome and genome analyses. Genome Res. 2009;19(4):521–32.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Mumbach MR, Rubin AJ, Flynn RA, Dai C, Khavari PA, Greenleaf WJ, Chang HY. HiChIP: efficient and sensitive analysis of protein-directed genome architecture. Nat Methods. 2016;13(11):919.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Jin F, Li Y, Dixon JR, Selvaraj S, Ye Z, Lee AY, Yen C-A, Schmitt AD, Espinoza CA, Ren B. A high-resolution map of the three-dimensional chromatin interactome in human cells. Nature. 2013;503(7475):290.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Bonev B, Cohen NM, Szabo Q, Fritsch L, Papadopoulos GL, Lubling Y, Xu X, Lv X, Hugnot, JP, Tanay A, Cavalli G. Multiscale 3D genome rewiring during mouse neural development. Cell 2017;171(3):557–72.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Stadhouders R, Vidal E, Serra F, Di Stefano B, Le Dily F, Quilez J, Gomez A, Collombet S, Berenguer C, Cuartero Y, Hecht J. Transcription factors orchestrate dynamic interplay between genome topology and gene regulation during cell reprogramming. Nat Genet 2018;50(2):238.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Browne G, Cipollone R, Lena AM, Serra V, Zhou H, van Bokhoven H, Dotsch V, Merico D, Mantovani R, Terrinoni A, et al. Differential altered stability and transcriptional activity of DeltaNp63 mutants in distinct ectodermal dysplasias. J Cell Sci. 2011;124(Pt 13):2200–7.

    Article  CAS  PubMed  Google Scholar 

  39. Celli J, Duijf P, Hamel BC, Bamshad M, Kramer B, Smits AP, Newbury-Ecob R, Hennekam RC, Van Buggenhout G, van Haeringen A, et al. Heterozygous germline mutations in the p53 homolog p63 are the cause of EEC syndrome. Cell. 1999;99(2):143–53.

    Article  CAS  PubMed  Google Scholar 

  40. Hong J-W, Hendrix DA, Levine MS. Shadow enhancers as a source of evolutionary novelty. Science. 2008;321(5894):1314.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Barolo S. Shadow enhancers: frequently asked questions about distributed cis-regulatory information and enhancer redundancy. BioEssays. 2012;34(2):135–41.

    Article  CAS  PubMed  Google Scholar 

  42. van Bokhoven H, Hamel BC, Bamshad M, Sangiorgi E, Gurrieri F, Duijf PH, Vanmolkot KR, van Beusekom E, van Beersum SE, Celli J, et al. p63 Gene mutations in eec syndrome, limb-mammary syndrome, and isolated split hand-split foot malformation suggest a genotype-phenotype correlation. Am J Hum Genet. 2001;69(3):481–92.

    Article  PubMed  PubMed Central  Google Scholar 

  43. Rheinwald JG, Green H. Epidermal growth factor and the multiplication of cultured human epidermal keratinocytes. Nature. 1977;265(5593):421–4.

    Article  CAS  PubMed  Google Scholar 

  44. Rinne T, Clements SE, Lamme E, Duijf PH, Bolat E, Meijer R, Scheffer H, Rosser E, Tan TY, McGrath JA, et al. A novel translation re-initiation mechanism for the p63 gene revealed by amino-terminal truncating mutations in Rapp-Hodgkin/Hay-Wells-like syndromes. Hum Mol Genet. 2008;17(13):1968–77.

    Article  CAS  PubMed  Google Scholar 

  45. Liu NQ, Ter Huurne M, Nguyen LN, Peng T, Wang S-Y, Studd JB, Joshi O, Ongen H, Bramsen JB, Yan J. The non-coding variant rs1800734 enhances DCLK3 expression through long-range interaction and promotes colorectal cancer progression. Nat Commun. 2017;8:14418.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Novakovic B, Habibi E, Wang S-Y, Arts RJ, Davar R, Megchelenbrink W, Kim B, Kuznetsova T, Kox M, Zwaag J. β-glucan reverses the epigenetic state of LPS-induced immunological tolerance. Cell. 2016;167(5):1354–68.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M, Gingeras TR. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29(1):15–21.

    Article  CAS  PubMed  Google Scholar 

  48. Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010;11(10):R106.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Li H, Durbin R. Fast and accurate short read alignment with Burrows–Wheeler transform. Bioinformatics. 2009;25(14):1754–60.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Zhang Y, Liu T, Meyer CA, Eeckhoute J, Johnson DS, Bernstein BE, Nussbaum C, Myers RM, Brown M, Li W, et al. Model-based analysis of ChIP-Seq (MACS). Genome Biol. 2008;9(9):R137.

    Article  PubMed  PubMed Central  Google Scholar 

  51. Ramírez F, Dündar F, Diehl S, Grüning BA, Manke T. deepTools: a flexible platform for exploring deep-sequencing data. Nucleic Acids Res. 2014;42(W1):W187–91.

    Article  PubMed  PubMed Central  Google Scholar 

  52. Tripathi S, Pohl MO, Zhou Y, Rodriguez-Frandsen A, Wang G, Stein DA, Moulton HM, DeJesus P, Che J, Mulder LC. Meta-and orthogonal integration of influenza “OMICs” data defines a role for UBR4 in virus budding. Cell Host Microbe. 2015;18(6):723–35.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Wingett S, Ewels P, Furlan-Magaril M, Nagano T, Schoenfelder S, Fraser P et al. HiCUP: pipeline for mapping and processing Hi-C data. F1000Res 2015; 4:1310.

    Article  PubMed  PubMed Central  Google Scholar 

  54. Cairns J, Freire-Pritchett P, Wingett SW, Várnai C, Dimond A, Plagnol V, Zerbino D, Schoenfelder S, Javierre B-M, Osborne C. CHiCAGO: robust detection of DNA looping interactions in Capture Hi-C data. Genome Biol. 2016;17(1):127.

    Article  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We thank Eva Janssen-Megens, Siebe van Genesen and Rita Bylsma for operating the Illumina analyzer. We thank the ENCODE Consortium for sharing their data.

Authors' contributions

JQ, GY and HZ conceived and designed the experiments. JQ, GY, HZ wrote and revised the manuscript. JQ performed the experiments. JQ, GY and HZ analyzed the data. All authors read and approved the final manuscript.

Funding

This research was supported by Netherlands Organization for Scientific Research (NWO/ALW/MEERVOUD/836.12.010, HZ), Radboud University fellowship (HZ) and Chinese Scholarship Council grant 201406330059 (JQ).

Data availability

To review our complete dataset GEO accession GSE123711 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE123711). All data supporting the findings of the study and in-house codes are available on request.

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Author information

Author notes

  1. Guoqiang Yi

    Present address: Center for Animal Genomics, Agricultural Genome Institute at Shenzhen, Chinese Academy of Agricultural Sciences, Shenzhen, 518124, China

  2. Jieqiong Qu and Guoqiang Yi contributed to the manuscript equally

Authors and Affiliations

  1. Department of Molecular Developmental Biology, Faculty of Science, Radboud Institute for Molecular Life Sciences, Radboud University, Nijmegen, The Netherlands

    Jieqiong Qu & Huiqing Zhou

  2. Department of Molecular Biology, Faculty of Science, Radboud Institute for Molecular Life Sciences, Radboud University, Nijmegen, The Netherlands

    Guoqiang Yi

  3. Department of Human Genetics, Radboud University Medical Centre, Nijmegen, The Netherlands

    Huiqing Zhou

Authors
  1. Jieqiong Qu
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Guoqiang Yi
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Huiqing Zhou
    View author publications

    You can also search for this author in PubMed Google Scholar

Corresponding author

Correspondence to Huiqing Zhou.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1: Figure S1.

Few overlap between p63 and CTCF binding sites. A Overlapping of p63 and CTCF binding sites from ChIP-seq data. B Heatmaps showing the CTCF binding signal at all 45,350 p63 binding sites (left) and the p63 binding signal at all 12,116 CTCF binding sites (right). C Bandplots showing the p63 binding signal (left) and CTCF binding signal (right) at the 435 co-binding sites in both control and p63 mutant keratinocytes.

Additional file 2: Figure S2.

Deregulated nucleosome organization at p63-bound and CTCF-bound Ctr-OCRs in p63 mutant keratinocytes. A Reads density of ATAC-seq at p63-bound Ctr-OCRs with different insert sizes in both control and p63 mutant keratinocytes. B Reads density of ATAC-seq at CTCF-bound Ctr-OCRs with different insert sizes in both control and p63 mutant keratinocytes.

Additional file 3: Figure S3.

CTCF bound preferentially to Ctr-OCRs marked by H3K4me3. Ctr-OCRs are grouped into two categories according to H3K4me3 signal, Ctr-OCRs with H3K4me3 (purple) and Ctr-OCRs without H3K4me3 (gray). CTCF ChIP-seq showed higher signal in Ctr-OCRs with H3K4me3, while p63 ChIP-seq showed higher signal in Ctr-OCRs without H3K4me3 in control keratinocytes. ***, p < 0.001, Wilcoxon signed-rank test.

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Qu, J., Yi, G. & Zhou, H. p63 cooperates with CTCF to modulate chromatin architecture in skin keratinocytes. Epigenetics & Chromatin 12, 31 (2019). https://doi.org/10.1186/s13072-019-0280-y

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13072-019-0280-y

Epigenetics & Chromatin

ISSN: 1756-8935