Search
Close this search box.

The transcriptional regulatory network modulating human trophoblast stem cells to extravillous trophoblast differentiation – Nature Communications

Transcriptome and enhancer dynamics during EVT differentiation

As little is known about EVT differentiation, we first carried out RNA-sequencing (RNA-seq) during time course differentiation of human TSCs to EVTs (Fig. 1a). Upon differentiation, epithelial to mesenchyme-like morphological changes and a corresponding increase in protein levels of HLA-G, an EVT marker gene, and a decrease in the expression of TP63, a TSC marker gene, were observed, indicating proper EVT differentiation (Fig. 1a and Supplementary Fig. 1a). The enhanced surface expression of HLA-G (>94% of HLA-G positive cell) further reinforced this observation (Supplementary Fig. 1b). The activation of EVT marker genes (HLA-G and MMP2) and downregulation of TSC markers (TP63 and TEAD4) were also observed in RT-qPCR and RNA-seq results (Fig. 1b and Supplementary Data 1). In addition, Gene Set Enrichment Analysis (GSEA) using the EVT and cytotrophoblast (CT, in vivo counterpart of TSCs) marker genes defined from the Molecular Signatures Database17,18 and a scRNA-seq study of the human first-trimester placenta19 further supported the proper upregulation of EVT-active genes along with the downregulation of CT genes, substantiating the proper EVT differentiation (Supplementary Fig. 1c). Principle component analysis (PCA) confirmed the gradual alteration of the transcriptome during EVT differentiation, which was distinct from the transcriptome of fully differentiated STs (Fig. 1c). The number of differentially expressed genes (DEGs; absolute log2-fold change >1 and p < 0.05) increased as differentiation progressed (Supplementary Fig. 1d), and we obtained a total of 5679 DEGs (corresponding to at least one of the time-points).

Fig. 1: Dynamic changes in transcriptome and enhancer landscape during EVT differentiation.
figure 1

a Representative brightfield images of TSCs and differentiating cells to EVTs on days 1 (EVT D1), 2 (EVT D2), 3 (EVT D3), 4 (EVT D4), 5 (EVT D5), and 8 (EVT D8). Scale bar: 100 µm. Two independent repeats were carried out, resulting in similar results. Cells used for transcriptome profiling (top) and enhancer mapping (bottom) are labeled. b Relative mRNA expression levels of EVT marker genes (HLA-G and MMP2) and TSC marker genes (TP63 and TEAD4) in TSCs and EVT differentiating cells. The fold change relative to TSCs is plotted. Error bars: mean ± SD (n = 3, independent repeats). c Principle component analysis illustrating the transcriptome profiles of TSCs, EVT differentiating cells, and STs. d Heatmaps showing Dirichlet Process Gaussian process (DPGP) mixture model clustering20 of genes that exhibited differential expression during at least one time point compared to TSCs (left), and H3K27ac signals assigned to classes 1–4 genes (right). The values represent the log2-transformed counts and H3K27ac signals, respectively, compared to the average. e GO analysis of classes 1–4 genes was conducted and p-values were calculated in Metascape89. f Sample-to-sample distances heatmap showing an overview of similarities and dissimilarities in H3K27ac signals among TSCs, EVT D3, and EVT D8.

To functionally characterize the DEGs, we performed a Dirichlet Process Gaussian process (DPGP) mixture model clustering20, and identified four distinct expression patterns (Fig. 1d and Supplementary Data 2). The genes belonging to class 1 showed high expression in TSCs and a gradual decrease in expression throughout differentiation (TSC-active). In contrast, class 2 genes were inactive in TSCs but continuously activated upon differentiation (EVT-active). Class 3 genes showed dynamic expression with strongest activity at the early-stage of EVT differentiation (at days 2 and 3) followed by rapid repression in mature EVTs (early-stage-active). Finally, class 4 genes were active in both TSCs and mature EVTs, but downregulated in between these two states.

To map corresponding changes in enhancer landscape, we performed chromatin immunoprecipitation coupled with next-generation sequencing (ChIP-seq) of H3K27ac21 during EVT differentiation at day 0 (TSC), early-stage differentiation at day 3 (EVT D3), and mature EVT at day 8 (EVT D8) (Fig. 1a). We identified 46,591, 44,056, and 42,666 H3K27ac-positive sites in TSC, EVT D3, and EVT D8, respectively (Supplementary Fig. 1e), confirming substantial and dynamic changes in enhancer usage during EVT differentiation (Fig. 1f). A large portion of these sites localizes within intergenic or introns as expected (Supplementary Fig. 1e). We further validated the enhancers identified by H3K27ac signals using other well-established enhancer markers, including P300, MED1, and H3K4me122. Supplementary Fig. 1f illustrates the overall similarities in the enhancer landscape. Specifically, we observed co-occupancy of P300 and MED1 at the enhancer loci, along with the enriched signals of H3K27ac and H3K4me1 close to these binding sites.

We confirmed that changes in enhancer landscape and gene expression are overall positively correlated, as observed in other contexts (Fig. 1d and Supplementary Fig. 1g)23,24,25,26. These results align with the notion that cell-type-specific enhancers determine cell-type-specific gene expression programs27,28,29. Validating our results, previously known TSC markers such as TP63, ELF5, and TEAD4 belonged to class 1, while the EVT markers, such as HLA-G and MMP2, and a known EVT regulator, ASCL216, were included in class 2 (Supplementary Fig. 1h and Supplementary Data 2). Furthermore, gene ontology (GO) analysis identified terms consistent with the known functions of each cell type. Class 1 genes were enriched in terms including cell division, cell cycle, and replication, consistent with rapidly proliferative TSCs11,12. Classes 2 and 4 genes were enriched in terms associated with EVT functions, such as invasive characteristics30 and interaction with immune cells, which is known to occur during spiral artery remodeling31. Class 3 genes were enriched in cellular response to stress and intrinsic apoptotic signaling pathways, consistent with the prior reports describing that differentiation occurs along with cellular stress and apoptosis32,33,34 (Fig. 1e).

Two classes of EVT regulators

As master TFs are often associated with clusters of enhancers known as super-enhancers (SEs)35, we determined SEs to identify candidate key regulators of EVT differentiation. We applied the ROSE algorithm35,36 to the H3K27ac ChIP-seq data and mapped SE-associated genes in TSC, EVT D3, and EVT D8 cells (Supplementary Fig. 2a and Supplementary Data 3). As anticipated, the previously known TSC-associated TFs TEAD412, YAP137, and MSX214 were associated with SEs in TSCs. A recently reported EVT regulator ASCL216 was associated with the SEs defined in EVT D8 cells (Supplementary Fig. 2a and Supplementary Data 3), confirming the validity of our approach.

We first focused on the genes belonging to the classes 2 and 4, which are highly active in fully differentiated EVTs (Fig. 1d), as these genes were also enriched in EVT-related GO terms (Fig. 1e). Among 394 class 2 genes, 18 were TFs38, and 7 TFs (DLX5, DLX6, ASCL2, ZNF439, IRF7, SNAI1, and VAX2) were associated with the SEs defined in EVT D8 (Fig. 2a–c). Among them, we selected 4 TFs (DLX5, DLX6, ASCL2, and ZNF439, hereafter, late-stage TFs) for further validation. We also included NRIP1 from class 4 as a late-stage TF candidate due to its association with multiple invasive cancers39,40 and SEs defined in EVT D8. NRIP1 also highly expressed in mature EVTs, similar to other late-stage TFs (Fig. 2b, c).

Fig. 2: Expression and enhancer association define two classes of EVT regulators.
figure 2

a Number of genes (blue), TFs (green), and SE-associated TFs (red) in class 2 and 3 genes. b Gene track view of H3K27ac signals near the candidate regulators of EVT differentiation, with red bars indicating SEs. c Expression patterns of early- and late-stage TF candidates during EVT differentiation. Values represent log2-transformed counts compared to the average. d Enriched motifs within the top 30% of H3K27ac peaks in EVT on differentiation day 3 were identified using the findMotifsGenome.pl module under HOMER (v4.11). p-values were calculated based on the default settings of the module, with a region size of 500 bp. e Western blot analysis demonstrating the expression patterns of EVT regulator candidates in CT27 TSC line. HLA-G expression was measured to validate EVT differentiation. ACTB was used as a loading control. *ZNF439 (lower band). Independent repeats, utilizing other TSC lines, were conducted as presented in Supplementary Fig. 2d. f Normalized transcripts per million value (nTPM) of EVT regulator candidates in scRNA-seq data of human first-trimester placenta (retrieved from the Human Protein Atlas19,102).

We surmised that TFs belonging to class 3 may also play roles in EVT differentiation, perhaps acting as early-stage regulators as they showed a unique expression pattern, peaking at 2–3 days of differentiation, followed by rapid downregulation (Fig. 1d). Since there were many TFs-associated with SEs in this group (39 TFs, Fig. 2a), we additionally conducted a motif analysis of the top 30% enhancer peaks defined in EVT D3 cells. Interestingly, TFAP2C, whose motif is similar to that of TFAP2A, was the only SE-associated TF that overlapped with the top 5 motifs identified (Fig. 2d). In mice, Tfap2c activates TSC self-renewal genes41 and has also been implicated in the trans-differentiation of mouse PSCs and fibroblasts into TSLCs42,43. Similarly, recent studies using human TSCs reported that TFAP2C is required for self-renewal of TSCs14,44. Intriguingly, TFAP2C has also been implicated in human epidermal lineage commitment as an initiation-stage regulator of keratinocyte specification45. Based on our data and the known roles of TFAP2C during cell fate commitment, we selected TFAP2C as an early-stage regulator candidate controlling EVT differentiation.

To validate the expression patterns of the identified candidates for EVT regulators, we conducted RT-qPCR and Western blot analyses in TSCs and differentiating cells towards EVT (CT27 TSC line). These experiments confirmed that the candidate EVT regulators displayed similar expression patterns to those observed in the RNA-seq results (Fig. 2e and Supplementary Fig. 2b). Notably, NRIP1 exhibited comparable expression levels in TSCs and EVTs on differentiation day 8 in the RNA-seq data. However, the expression pattern of NRIP1 during EVT differentiation, as observed in the RT-qPCR and Western blot results, resembled that of other late-stage TFs (Fig. 2e and Supplementary Fig. 2b, d). We further confirmed the consistent morphological changes and expression patterns of the candidate EVT regulators during EVT differentiation in other TSC lines (CT29 and CT30; Supplementary Fig. 2b–d).

To further determine the in vivo relevance of the candidate TFs, we verified the expression of selected TFs in human placentas. Data from single-cell RNA-seq (scRNA-seq) studies of the human first-trimester placentas19,46 confirmed that the late-stage TFs, DLX5, DLX6, ASCL2, and NRIP1 are also highly expressed in EVTs (Fig. 2f and Supplementary Fig. 2e, f). As expected, TFAP2C expression was strong in CT and lower in EVT populations (Fig. 2f and Supplementary Fig. 2e, f). Overall, existing data validate our approach and classification scheme. Our transcriptome and enhancer analysis therefore identified a set of regulators expressed both early (TFAP2C) and late (DLX5, DLX6, ASCL2, ZNF439, and NRIP1) stages that were candidates for regulating the EVT differentiation process.

Early and late-stage TFs are essential for EVT differentiation

To determine whether these candidate TFs are required for proper EVT differentiation, we conducted short hairpin RNA (shRNA)-mediated knockdown (KD) in TSCs and subsequently differentiated them to EVTs (Fig. 3a and Supplementary Fig. 3a). While control TSCs showed elongated and spindle-like cell morphology upon differentiation, KD of both early- and late-stage TFs led to an abnormal cellular morphology (Supplementary Fig. 3b) and the failure to induce EVT markers, such as HLA-G and MMP2 (Fig. 3a and Supplementary Fig. 3c). Notably, the impaired induction of the EVT markers was also observed in CT29 and CT30 TSC lines upon KD of the candidate TFs (Supplementary Fig. 3a, c). Since one critical characteristic of EVTs is their invasion capacity, we examined the invasiveness of the KD cells using invasion chambers (Supplementary Fig. 3d)47. Compared to control EVTs, all KD cells showed significantly reduced invasion ability, indicating that both early- and late-stage TFs are essential for normal EVT differentiation (Fig. 3b).

Fig. 3: Both early and late-stage TFs are required for EVT differentiation.
figure 3

a Western blot analysis showing KD efficiency and expression level of EVT marker genes in EVTs following KD of each candidate EVT regulator. EVTs differentiated from TSCs infected with lentivirus expressing non-targeting sequence were used as control. Two independent repeats were carried out, resulting in similar results. b Invasion ability of EVTs after KD of each candidate EVT regulator, compared to control cells. Scale bar: 100 µm. The invasion area was quantified using Image J software. Measurements were obtained from images taken in three independent biological replicates. Error bars: mean ± SD (n = 3, independent repeats). c Heatmap showing relative expression of class 1–4 genes in EVT factor KD cells compared to non-targeting control. Color represents log2 mRNA expression relative to control KD. d GSEA using EVT marker genes defined in scRNA-seq data of the human first-trimester placenta19,102 and fetal placenta EVT gene set obtained from the Molecular Signatures Database, comparing EVT factor KD cells with control. NES, normalized enrichment score. FDR, false discovery rate. e GO analysis was performed on genes exhibiting lower expression in EVTs after KD of each EVT factor compared to the control. The analysis was conducted using DAVID88, and p-values were calculated. Color representation indicates significance, with -log(p-value). f Principal component analysis of EVT factor KD and control cell transcriptome. g GO analysis was performed on genes exhibiting higher expression in EVTs after KD of each EVT factor compared to the control. The analysis was conducted using DAVID88, and p-values were calculated. Color representation indicates significance, with -log(p-value).

To further explore the importance of the early- and late-stage TFs in controlling EVT-specific gene expression programs, we profiled the transcriptome of each TF KD cell line upon EVT differentiation by bulk RNA-seq. We identified 1534-2072 DEGs in KD cells compared to the control cells (Supplementary Fig. 3e, f and Supplementary Data 4). Firstly, we monitored the expression of the class 1–4 genes we defined earlier (Fig. 1d). As shown in Fig. 3c, KD of these TFs led to defective induction of genes active in EVTs (classes 2 and 4) compared to control. To further validate the results, we performed GSEA with the EVT-active gene set defined from the scRNA-seq of the human placenta by another group19, confirming the lower expression of the EVT signature genes in each TF KD line. The GSEA using another EVT signature gene set obtained from the Molecular Signatures Database17,48 also generated similar results (Fig. 3d). Consistently, metallopeptidases, such as MMP249,50 and ADAM1951,52, as well as previously known regulators of invasion, such as CXCR453,54 and EGFR-AS155, were downregulated in the EVT factor KD cells. In addition, NOTUM, an extracellular Wnt deacylase56,57, was commonly downregulated in the EVT factor KD cells, suggesting impaired EVT differentiation, as Wnt inhibition is necessary for EVT lineage formation11,58. We also observed a decrease in the expression of genes associated with interstitial EVT (iEVT) differentiation, such as ITGA1, PLAC8, and SREPINE259, in the EVT factor KD cells. To gain further insight into the function of each TF, we conducted GO analysis with the genes showing lower expression in the KD cells compared to the control cells (Supplementary Data 4). We found that the genes showing lower expression in both early- and late-stage TF KD cells are enriched in terms previously implicated in EVT characteristics and functions, such as extracellular matrix organization, cell migration, extracellular matrix disassembly, Jak-STAT and PI3K-Akt signaling pathways60, and response to hypoxia61 (Fig. 3e). We concluded that both TFAP2C and the late-stage TFs are essential for normal EVT function and induction of EVT-active genes.

TFAP2C KD leads to a unique expression pattern

Although the KD of all individual TFs triggers similar differentiation defects, we surmised that there might be functional differences between early- and late-stage TFs. As shown in Fig. 3c, we found that TFAP2C KD cells showed distinct gene expression patterns among the tested TFs. While KD of late-stage TFs showed improper downregulation of class 1 (TSC-active) genes, overall expression of the class 1 genes in TFAP2C KD cells was even lower than that in the control cells. This implies that early- and late-stage TFs may have unique roles in regulating class 1 genes. Since TFAP2C is significantly expressed in self-renewing TSCs, TFAP2C may act as an activator of class 1 genes, and it may play important roles in not only EVT differentiation but also TSC self-renewal as recently suggested14,44. On the other hand, KD of the late-stage TFs resulted in failed downregulation of class 1 genes, implying that the repression of TSC-active genes does not occur efficiently when these factors are knocked down. This observation was further supported by the results of GSEA using the CT marker genes defined from the placenta scRNA-seq study19 (Supplementary Fig. 3g). Similarly, several genes associated with cell cycle regulation, such as CDK1 and CCND1, along with the factors involved in cell proliferation, including MYBL262 and TOP2A (DNA topoisomerase 2-alpha, which plays a role in DNA replication and cell division), were significantly upregulated in the late-stage TF KD cells. Moreover, the late-stage TF KD cells exhibited sustained expression of the genes associated with epithelial cells, including EPCAM63, CDH164, and TJP1 (tight junction protein 1). Additionally, the key genes involved in Wnt signaling, such as AXIN2, FZD5, LRP5, and TCF7L165,66,67, displayed increased expression in the late-stage TF KD cells, implying that the KD cells still retained a TSC-like gene expression program (Supplementary Data 4). Thus, late-stage TFs might be involved in the repression of TSC-active genes directly or indirectly during differentiation. We also examined the changes in ST marker genes19 in the KD cells. Consistent with the previous report16, we observed increased expression of ST lineage genes in ASCL2 KD cells compared to control cells (Supplementary Fig. 3h). Likewise, TFAP2C KD, DLX6 KD, and NRIP1 KD cells exhibited higher expression of ST genes compared to the individual controls. Despite the increased expression of the ST gene set in TFAP2C KD and some late-stage TF-KD cells, the overall differences in the expression between TFAP2C KD and late-stage TF KD cells were confirmed by PCA (Fig. 3f). The late-stage TF KD cells were clustered together but separated from the TFAP2C KD or control cell clusters.

To further understand the TF-specific differences in gene regulation, we conducted GO analysis with the genes showing higher expression in the KD cells than in control cells (Supplementary Data 4). Consistent with Fig. 3c, the genes that failed to be repressed upon KD of late-stage TFs were enriched in cell proliferation and Hippo signaling pathway, the terms associated with proliferation of TSC and inhibition of differentiation12,37. However, this was not the case in TFAP2C KD cells (Fig. 3g). We concluded that while both TFAP2C and late-stage TFs are commonly required for proper EVT differentiation and the activation of EVT-active genes, they also play unique roles during EVT differentiation, likely in a differentiation stage-specific manner.

Critical role of early TFAP2C expression in EVT differentiation

TFAP2C shows a unique expression pattern during EVT differentiation, peaking early and then decreasing (Fig. 2c). To gain further insight into the role of this TF, and to test the significance of its expression dynamics, we employed both loss- and gain-of-function studies (Fig. 4a). First, to test if TFAP2C plays distinct roles depending on the stage of differentiation, we performed KD of TFAP2C at different time points during differentiation. Similar to the previous results, KD of TFAP2C at the early stage (KD at day −1, KD-Early) resulted in impaired EVT differentiation evidenced by defects in morphological changes, EVT marker gene activation, and invasion ability. Intriguingly, we did not detect such defects upon KD of TFAP2C beginning at day 3 (KD-Mid) or 5 (KD-Late) (Fig. 4b–d and Supplementary Fig. 4a), which was verified with multiple EVT marker genes (Supplementary Fig. 4b).

Fig. 4: Stage-specific regulation of TFAP2C during EVT differentiation.
figure 4

a Experimental design for time-course TFAP2C KD (top) and OE (bottom). b Relative expression of TFAP2C and EVT marker genes in EVTs after KD of TFAP2C at different time points. E: KD day −1–8, M: KD day 3–8, L: KD day 5–8. Error bars: mean ± SD (n = 3, independent repeats). Significance by two-sided Student’s t-test. c Invasion ability of EVTs after KD of TFAP2C at different time points and individual control cells. Early: KD day −1–8, Mid: KD day 3–8, Late: KD day 5–8. Scale bar: 100 µm. d Immunofluorescence of EVT marker genes in EVTs after TFAP2C KD at different time points. Scale bar: 100 µm. Two independent repeats were carried out, resulting in similar results. e Brightfield images of EVTs overexpressing TFAP2C (Dox-treated day 2–8) and control. Scale bars: 100 µm. Two independent repeats were carried out, resulting in similar results. f Relative expression of EVT markers in TFAP2C OE compared to each control. Dox was treated on day 2, and cells were collected on days 3, 5, and 8. Error bars: mean ± SD (n = 3, independent repeats). Significance by two-sided Student’s t-test. g Immunofluorescence of HLA-G and TP63 in TFAP2C OE (day 2–8) and control. Scale bar: 100 µm. Two independent repeats were carried out, resulting in similar results. h Box plot illustrates class 1–4 gene expression in TFAP2C OE (day 2–8) compared to control (***p < 2.2e-16, two-tailed Wilcoxon rank-sum test). Each group has two independent RNA-seq results. The median (black bar) and box (25th to 75th percentiles) represent data distribution, while whiskers denote the minimum to maximum values. i, GSEA of TFAP2C OE and control using the EVT and CT marker gene sets defined in scRNA-seq data of human first-trimester placenta19,102. NES normalized enrichment score, FDR false discovery rate. j Correlation heatmap showing comparisons of TFAP2C OE with EVTs after KD of individual EVT regulators. Spearman correlations were calculated using log2-fold changes of the class 1–4 genes.

As TFAP2C is essential for TSC maintenance and the KD of TFAP2C in TSCs induces apoptosis after four days68, it is possible that TFAP2C KD in TSCs could impact downstream cellular physiology, thereby impeding EVT formation. To differentiate between the effects of TFAP2C KD in self-renewing TSCs and during differentiation, we performed additional KD experiments at multiple time points during the early stage of EVT differentiation. As shown in Supplementary Fig. 4c, the cells treated with lentivirus targeting TFAP2C on days −1, 0, and 1 exhibited similar impairments in the induction of EVT marker genes. However, TFAP2C KD on days 2 and 3 did not lead to defects in EVT marker gene induction. Additionally, we observed that apoptosis was not significantly induced in TFAP2C KD cells after 18 h of lentivirus infection, which aligns with the timing of EVT differentiation initiation (Supplementary Fig. 4d). These findings indicate that the depletion of TFAP2C in TSCs may not be the primary cause of defects in EVT differentiation. Instead, TFAP2C plays a crucial role during the early stage of EVT differentiation. In contrast, the depletion of late-stage TFs on days 3 and 5 hindered the induction of EVT marker genes, suggesting that these late-stage TFs are necessary for activating EVT-active genes during the late stage of differentiation (Supplementary Fig. 4e).

Next, we tested whether the downregulation of TFAP2C during the late stage of differentiation is also important for EVT differentiation. To do so, we maintained TFAP2C expression beyond the early stage of differentiation. We over-expressed TFAP2C via a doxycycline (Dox)-inducible system (see methods section and Supplementary Fig. 4f) on day 2 of differentiation (TFAP2C OE), when the endogenous TFAP2C level peaks and starts decreasing (Figs. 2c and 4a). Although the ectopic expression of TFAP2C was not sustained till day 8, we observed that ectopic TFAP2C induces morphological abnormalities and defects in TSC marker repression as well as EVT marker activation (Fig. 4e–g, and Supplementary Fig. 4g–j). Accordingly, transcriptomic analysis revealed that EVT-active genes (classes 2 and 4) are not adequately upregulated, while the expression of TSC-active genes (class 1) is not downregulated normally in the TFAP2C OE cells (Fig. 4h). This observation is further substantiated by GSEA results employing previously defined EVT and CT gene sets19 (Fig. 4i). Interestingly, the pattern of expression of classes 1–4 genes in the TFAP2C OE cells was most similar to the pattern in cells lacking late-stage TFs cells (Fig. 4j), suggesting that prolonged expression of TFAP2C has a similar effect to loss of late-stage TFs. In conclusion, our data reveal that an optimal level and timing of TFAP2C expression are critical for proper EVT differentiation.

TFAP2C primes late-stage TFs and EVT-active genes

Our data thus far show that, despite its requirement for the activation of EVT-active genes, TFAP2C is primarily essential during the early stage of differentiation when EVT-specific genes are not yet transcriptionally active. Moreover, TFAP2C needs to be downregulated during the late stage of differentiation for proper induction of EVT-active genes. To investigate the mechanistic basis of these effects, we examined the relationship of TFAP2C to DLX6, one of the late-stage TFs, using time-course ChIP-seq (Fig. 5a). Similar to its expression pattern during EVT differentiation (Fig. 2c and Supplementary Fig. 5a, b), TFAP2C bound to its targets in self-renewing TSCs and cells in early-stage day 2, followed by decreased occupancy as differentiation progressed. On the other hand, we detected DLX6 target occupancy starting at day 5 of differentiation (Fig. 5b and Supplementary Fig. 5c, d). As expected, the TFAP2C motif was the most enriched motif in the TFAP2C binding loci (Supplementary Fig. 5e), and the DLX3 motif, which is similar to the DLX5 or DLX6 motif (Supplementary Fig. 5f), was enriched in the DLX6 target loci (Fig. 5c). Interestingly, the TFAP2C motif was also significantly enriched in DLX6 target loci, suggesting common targets for TFAPC2 and the late-stage TF DLX6 (Fig. 5c).

Fig. 5: TFAP2C primes late-stage TFs and EVT-active genes at the early-stage of differentiation.
figure 5

a Experimental design for time-course ChIP-seq of TFAP2C and DLX6 during EVT differentiation. b Number of TFAP2C and DLX6 ChIP-seq peaks in TSCs and differentiating cells to EVT on days 2, 5, and 8. c Enriched motifs in DLX6 ChIP-seq peaks in EVT day 5 (left) and 8 (right). Motifs were identified using the findMotifsGenome.pl module under HOMER (v4.11), and p-values were calculated using the module’s default settings. d Gene track view of TFAP2C (blue) and DLX6 (red) signals in TSCs, EVT D2, EVT D5, and EVT D8 at the locus near HLA-G. Heatmaps of (e) TFAP2C and (f) DLX6 signals on TFAP2C in TSCs unique target loci (G1), TFAP2C/DLX6-common target loci (G2), and DLX6 in EVT D8 unique target loci (G3) in TSCs, EVT D2, EVT D5, and EVT D8. PC, peak center. Heatmaps of TFAP2C in TSCs and EVT D2, and DLX6 in EVT D8 near the genes in (g) class 2 and (h) class 4. i Occupancy score of TFAP2C (left) and DLX6 (right) on DLX6 in EVT D8 unique target loci (G3) in TSCs, EVT D2, EVT D5, and EVT D8.

The overlapping targets between TFAP2C and DLX6 (Fig. 5d and Supplementary Fig. 5g) suggested a potential regulatory relationship between these factors, and perhaps between TFAP2C and other late-acting genes. To further explore this relationship, we conducted a differential binding analysis comparing target loci of TFAP2C in TSCs and DLX6 in EVTs. This analysis allowed us to identify three distinct groups of loci: unique target loci of TFAP2C in TSCs (Group 1, G1), TFAP2C-DLX6 common target loci (Group 2, G2), and unique target loci of DLX6 in EVTs (Group 3, G3; Supplementary Fig. 5h). We then plotted the ChIP-seq peaks of TFAP2C and DLX6 according to the differentiation stage and target groups (Fig. 5e, f). Each row in the plot represents the same genomic locus, allowing for observing the changes in TFAP2C and DLX6 binding patterns at different time points. As anticipated, G1 and G3 peaks are associated with TSC- and EVT-active genes, respectively (Supplementary Data 5). Intriguingly, many EVT-active genes were also found in G2 loci (Supplementary Data 5), suggesting that, at different differentiation stages, both TFAP2C and DLX6 occupy regulatory elements of the common EVT-active genes. We observed that matrix metalloproteinases, specifically MMP2 and MMP9, which are essential for invasion into the extracellular matrix, are associated with G2 loci. Additionally, we identified master regulators of epithelial-to-mesenchymal transition (EMT), namely SNAI1 and TWIST169, as targets of both TFAP2C and DLX6. This suggests that TFAP2C and DLX6 play crucial roles in controlling the genes responsible for the invasive function of EVTs. Notably, we also noticed that the angiogenic factor VEGF and its receptors KDR and FLT170 are occupied by both TFAP2C and DLX6, suggesting that both early- and late-stage TFs are involved in the activation of genes associated with spiral artery remodeling, a critical function of EVT. Focusing on genes in classes 2 and 4, as defined in Fig. 1d, we also confirmed TFAP2C occupancy at regulatory loci near the EVT-active genes in both TSC and EVT D2. These loci are later occupied by DLX6 in EVT D8 (Fig. 5g, h). Interestingly, we also observed that the binding of TFAP2C at G3 loci increases during early differentiation (from day 0 to 2) (Fig. 5e, i). While some EVT-active genes are pre-occupied by TFAP2C in TSCs (G2), others show increased occupancy of TFAP2C in EVT D2 (G3). TFAP2C’s occupancy on these EVT-active genes during the late stage of differentiation was not obviously detectible as the level of TFAP2C greatly diminished. Since EVT-active genes are not significantly transcribed in TSCs and EVT D2 cells, our results suggest that TFAP2C primes EVT-active genes during the early stage of differentiation.

Given that the EVT factor KD cells exhibited similar impairments in the induction of EVT-active genes, we further investigated the relationships between TFAP2C and the late-stage TFs by conducting an integrative analysis of TFAP2C/DLX6 targets and the RNA-seq data from the individual EVT TFs KD cells. As shown in Supplementary Fig. 5i, approximately 75–84% of the DEGs upon each KD were directly bound by TFAP2C, DLX6, or both. Interestingly, around 49% and 60% of the upregulated and downregulated genes in each KD were common targets of TFAP2C and DLX6 (G2), suggesting functional similarities among the late-stage TFs and potential collaborative works between other the late-stage TFs and TFAP2C.

Potential pioneering activity of TFAP2C on EVT-active genes

The pre-occupancy of TFAP2C for the future activation of EVT-active genes is reminiscent of pioneer factor activity71,72 (Fig. 5d, g, h). Indeed, the patterns of TFAP2C and DLX6 binding to target loci near late-stage TFs and EVT markers (HLA-G and MMP2) are consistent with this role (Fig. 6a and Supplementary Fig. 6a, b). Our data show that TFAP2C occupies cis-regulatory elements of the late-stage TFs and EVT markers during the early stage of differentiation, and DLX6 binds to the loci initially primed by TFAP2C during the late stage (Fig. 6a and Supplementary Fig. 6a, b). To test whether TFAP2C functions as a pioneer factor, we assessed TFAP2C target loci in the context of the genome-wide chromatin accessibility landscape monitored by ATAC-seq. We also mapped H3K4me3 signatures, a hallmark of active promoters, in TSCs and EVTs, to monitor the activity of the TFAP2C target genes. We found that TFAP2C occupies both open and closed chromatin regions near late-stage TFs and EVT markers (Fig. 6a and Supplementary Fig. 6a, b). Notably, several genomic loci bound by TFAP2C were near-closed chromatin in TSCs but showed significantly increased accessibility in EVTs. These TFAP2C-primed regions also exhibited increased H3K27ac signals and increased transcriptional activity of the associated genes in mature EVTs.

Fig. 6: TFAP2C acts as a pioneer factor for late-stage TFs and EVT-active genes.
figure 6

a Gene track view displaying TFAP2C and DLX6 signals in TSCs and differentiating cells to EVT on days 2, 5, and 8, along with ATAC-seq, H3K27ac, H3K4me3 signals in TSCs and EVTs. mRNA expression levels are shown near EVT factors. b Heatmaps showing TFAP2C ChIP-seq unique, common, and ATAC-seq unique loci in TSCs. TFAP2C ChIP-seq unique loci are indicated by a black bar. PC peak center. c Enriched GO terms of the top 3000 TFAP2C unique, common, or ATAC-seq unique peaks were determined using GREAT software100. p-values were calculated within the software. d ATAC-qPCR results validating TFAP2C binding to loci near EVT-active genes with low ATAC-seq signal. To account for primer variability, the Ct value was first normalized by subtracting the Ct value of sheared genomic DNA, resulting in a normalized Ct value. The deltaCt value was then calculated by subtracting the normalized Ct value of a genomic locus near an olfactory receptor gene (OR1A1), which exhibits very low ATAC-seq signal and low transcription levels in TSCs. The y-axis represents negative deltaCt values, indicating that an increase in value corresponds to increased openness. Error bars: mean ± SD (n = 3, independent repeats). e ATAC-seq signals at TFAP2C unique and common loci identified by comparing ChIP-seq results of TFAP2C in TSCs and DLX6 in EVTs. Box plots represent the median (black bar) and interquartile range (box boundaries), with whiskers from minimum to maximum values. Significance by two-tailed Wilcoxon rank-sum test. Analysis involved ATAC-seq signals on all TFAP2C unique or common loci in TSCs and EVTs. f Enhancer signals at TFAP2C unique loci in TSCs and EVTs. Box plots represent the median (black bar) and interquartile range (box boundaries), with whiskers from minimum to maximum values. Significance by two-tailed Wilcoxon rank-sum test. Statistical analysis involved signals from enhancer ChIP-seq results on all specified TFAP2C unique loci in both TSCs and EVTs.

To further investigate the potential pioneer factor activity of TFAP2C during EVT specification, we compared TFAP2C binding loci and ATAC-seq signals in TSCs using the MAnorm73. We observed that TFAP2C binding loci with low ATAC-seq signals (Fig. 6b, indicated by a black bar), which are primarily found in enhancer regions (Supplementary Fig. 6c). The genes associated with these loci are enriched in the terms, such as JAK-STAT and PI3K-Akt, previously implicated in EVT differentiation (Fig. 6c). The vasculature development-related terms were also consistent with the spiral artery remodeling function of EVTs (Fig. 6c). These results indicate that TFAP2C binding loci with low accessibility are associated with EVT-active genes. To validate the binding of TFAP2C on genomic loci with low ATAC-seq signals near EVT-active genes in TSCs, we conducted an ATAC-qPCR analysis. We examined the loci with varying degrees of openness in TSCs, including those with very weak ATAC-seq signals near olfactory receptor genes (inactive), weak ATAC-seq signals near EVT-active genes not active in TSCs (EVT-active), strong ATAC-seq signals near TSC marker genes (TSC-active), and ribosomal protein-coding genes (TSC/EVT-active) (Supplementary Fig. 6d). The ATAC-qPCR results in Fig. 6d confirmed that the genomic loci near EVT-active genes bound by TFAP2C in TSCs (red bars) exhibit overall lower openness compared to the TSC-active loci (blue bars) and TSC/EVT-active loci (green bars), while still being more open than closed loci near olfactory receptor genes (inactive, black bars).

Consistently, the ATAC-seq signals at the TFAP2C unique loci were significantly increased in EVTs, along with increased enhancer signals (Fig. 6e, f). Thus, TFAP2C target loci with low accessibility during the early EVT differentiation become more accessible as differentiation proceeds, further supporting the idea of pioneering activity for TFAP2C. As KD of TFAP2C at the early differentiation stage inhibits EVT differentiation (Fig. 3), TFAP2C-dependent priming of EVT-active genes appears essential for EVT differentiation. Together, these results support the model that TFAP2C exhibits potential pioneer factor activity for specific targets during EVT differentiation.

Late-stage TFs form a transcriptional regulatory network

We then sought to determine the mechanisms by which late-stage TFs regulate EVT differentiation. DLX6 occupancy on the regulatory elements of EVT-active genes suggests direct transcriptional regulation. To further investigate the function of late-stage TFs, we mapped their genomic targets by ChIP-seq (DLX5 and ASCL2) and bioChIP-seq (NRIP1)74 in mature EVT D8 cells (Supplementary Fig. 4f). We could not detect significant ZNF439 binding peaks although we detected both nuclear and cytoplasmic localization of ZNF439 in EVT D8 cells (Supplementary Fig. 7a). While the number of binding loci varied (3413–86,735) depending on the TF (Supplementary Fig. 7b), motif analysis revealed that DLX and ASCL2 motifs are significantly enriched in DLX5 and ASCL2 target loci, respectively (Supplementary Fig. 7c). Similar to the motif analysis results obtained for DLX6 (Fig. 5c), we observed significant enrichment of the TFAP2C motif in binding loci of both DLX5 and ASCL2, implying a similar target priming mechanism by TFAP2C. An NRIP1 binding motif has not been reported, but we found enriched motifs for SP and KLF family proteins and TSC-related factors, such as TEAD3 and FOSL2 within the NRIP1 target loci (Supplementary Fig. 7c). Thus, we were able to map target occupancy patterns of most of the late TFs and confirm the validity of our results by motif analysis.

As cell-type specific key TFs often auto-regulate their own gene expression29,75, we asked whether the late-stage TFs also bind to their own regulatory regions. As shown in Fig. 7a and Supplementary Fig. 7d, DLX5, DLX6, and ASCL2 (except for NRIP1) occupied regulatory regions of late-stage TFs and EVT markers, forming auto-regulatory loops. Notably, DLX6, DLX5, and ASCL2 bound the regulatory elements of TFAP2C, suggesting potential roles of late-stage TFs in the regulation of TFAP2C expression during late stage of differentiation (Figs. 6a and 7a).

Fig. 7: Late-stage TFs collaborate in regulating the EVT-specific gene expression program.
figure 7

a Gene track view of late-stage TF binding loci in mature EVTs (EVT D8), near EVT factor loci. b Peak scores of late-stage TFs on class 1–4 genes. Box plots represent the median (black bar) and interquartile range (box boundaries), with whiskers indicating minimum to maximum values. Statistical significance by two-tailed Wilcoxon rank-sum test (NS not significant). All peak scores identified from ChIP-seq results were employed for statistical analysis. c Heatmaps of class 1–4 gene expression in EVTs after KD of individual late-stage TFs. Genes were ordered based on relative gene expression levels in the KD cells compared to the control (from lowest to highest). Occupancy signals plotted as a moving average (window size: 50). Red bars indicate enriched signals. d Relative mRNA expression of EVT marker genes in TSCs overexpressing individual late-stage TFs, induced by administrating Dox for 6 days. Error bars: mean ± SD (n = 3, independent repeats). Statistical significance was determined by two-sided Student’s t-test. e Venn diagram illustrating overlap of top 2000 target genes of DLX5, DLX6, and ASCL2. f Model: TFAP2C primes late-stage TFs during early stage of differentiation, while late-stage TFs activate themselves and repress TFAP2C, forming a transcriptional regulatory circuit to ensure proper EVT differentiation.

To obtain insights into the function of the late-stage TFs on their target genes, we examined TF occupancy on the four gene classes we defined earlier (Fig. 1d). As NRIP1 showed unique global occupancy and target gene regulation patterns (Fig. 7a and Supplementary Fig. 7d, e), we focused on DLX5, DLX6, and ASCL2. As shown in Fig. 7b, DLX5, DLX6, and ASCL2 showed significantly higher occupancy signals on the EVT-active genes (classes 2 and 4) compared to genes in classes 1 and 3, consistent with a role in activation of class 2 and 4 genes. To test late-stage TFs’ role as a direct activator of EVT-active genes, we performed combined analysis of occupancy and changes in expression upon KD of the late-stage TFs on class 1–4 genes (Fig. 7c). We observed increased occupancy patterns for the late-stage TFs at genes downregulated upon KD of each TF (activated genes by the TFs). Moreover, a majority of affected genes were from classes 2 and 4 (Fig. 7c). To further confirm this, we generated individual inducible OE cell lines of the late-stage TFs (Supplementary Fig. 7f) and performed OE of the late-stage TFs in self-renewing TSCs (Supplementary Fig. 7g). We observed the activation of EVT markers HLA-G and MMP2 (Fig. 7d).

As key TFs co-regulate cell-type specific genes by forming transcriptional regulatory networks29,74,76,77,78, we monitored the top 2000 target genes of DLX5, DLX6, and ASCL2 and confirmed that they indeed share many common targets (Fig. 7e). GO analysis of the common target genes revealed their enrichment in EVT function-related terms, such as Rho GTPase cycle, regulation of cell/focal adhesion, and tube morphogenesis (Supplementary Fig. 7h). Similar to our analysis shown in Fig. 7c, the common target genes were also the most significantly affected genes (classes 2 and 4 genes) upon individual KD of all late-stage TFs including ZNF439 and NRIP1 (Supplementary Fig. 7i). Altogether, we reveal that the late-stage TFs we identified form a transcriptional regulatory network and positively control EVT-active genes during EVT specification.

Latest Intelligence