Abstract The vascular endothelium features unique molecular and functional properties across different vessel types, such as between arteries, veins and capillaries, as well as between different organs, such as the leaky sinusoidal endothelium of the liver versus the impermeable vessels of the brain. However, the transcriptional networks governing endothelial organ specialization remain unclear. Here we profile the accessible chromatin and transcriptional landscapes of the endothelium from the mouse liver, lung, heart, kidney, brain and retina, across developmental time, to identify potential transcriptional regulators of endothelial heterogeneity. We then determine which of these putative regulators are conserved in human brain endothelial cells, and using single-cell transcriptomic profiling, we define which regulatory networks are active during brain maturation. Finally, we show that the putative transcriptional regulators identified by these three approaches molecularly and functionally reprogram naive endothelial cells. Thus, this resource can be used to identify potential transcriptional regulators controlling the establishment and maintenance of organ-specific endothelial specialization. __________________________________________________________________ Despite a shared mesodermal origin and a host of common functions, endothelial cells (ECs), which form the inner lining of all blood vessels, show extensive functional heterogeneity and phenotypic specialization^[35]1–[36]4. For example, the hepatic sinusoidal capillary vessels of the liver feature large intercellular gaps (or fenestrae) between ECs and they lack an organized basement membrane, which allow maximal contact and exchange between blood and hepatocytes^[37]5. These fenestrae are also essential for receptor-mediated endocytosis of lipoproteins, enabling sinusoidal ECs to function as scavengers, eliminating soluble macromolecular waste. By contrast, the primary function of ECs within the kidney glomeruli is to filter fluids and solutes^[38]6. While glomerular capillary ECs also possess inter cellular fenestrae, these gaps are smaller than in their liver sinusoidal counterparts and also cover more of the cell surface^[39]7. Unlike sinusoidal ECs, glomerular ECs secrete and deposit glycocalyx, a thick cell surface layer of membrane-associated proteoglycans, glycolipids, glycosamines and associated plasma proteins, which forms another filtration barrier^[40]8. In addition to these interorgan-specific features, ECs also feature extensive intraorgan diversity, as the endothelium within the same organ is endowed with unique characteristics. This heterogeneity is exemplified by the unique morphological features and molecular signatures that distinguish the endothelium of arteries, arterioles, veins, venules and capillary vessels^[41]9. In addition, single-cell transcriptional profiling studies have uncovered additional distinct EC subpopulations, such as Car4 ECs in the lung^[42]10 and up to 24 EC subtypes in the kidney^[43]11. These diverse functions and morphologic features demonstrate the inherent phenotypic plasticity of the endothelium and suggest that cellular heterogeneity is a core property that enables ECs to fulfill their sometimes quite disparate tasks in very different environments (the adaptations required to thrive in the harsh environs of the kidney are distinct from those needed in capillaries located proximal to alveoli within the oxygen-rich environment of the lung). However, all ECs share a common mesodermal origin and the same genome. Thus, their functional diversification probably derives from how that genome is activated via chromatin accessibility and/or epigenetic regulation^[44]12,[45]13. Enhancers, non-coding regions of the genome that modify transcriptional output, are central nodes in transcriptional networks, integrating multiple upstream signals into unified outputs that act to regulate promoter activity and ultimately induce changes in gene expression^[46]14. Transcriptional regulatory elements are difficult to predict a priori due to their undefined sequence or location (with respect to their target genes). However, while several techniques have emerged to map enhancers (for example, H3K27ac (histone 3, acetylated lysine 27) ChIP-seq, DNase hypersensitivity mapping, Ep300 chromatin occupancy), they each come with significant caveats. In addition, purifying endothelium from different organs is not trivial, and complicated fluorescence activated cell sorting (FACS) procedures probably introduce artifacts given the extensive delay between tissue collection and EC isolation. Furthermore, the amount of input material required can be daunting if the lineage of interest comprises a small fraction of the cells in a tissue (for example, the adult mouse retina, which only contains approximately 5,000 ECs). By combining Cre-dependent expression of a genetically encoded, fluorescent and epitope-tagged nuclear membrane protein (Sun1-2xsfGFP-6xMyc)^[47]15 with an endothelial-specific CreER driver line^[48]16, we rapidly isolated endothelial nuclei by affinity pulldown via isolation of nuclei tagged in specific cell types (INTACT)^[49]17. By processing a fraction of endothelial nuclei for the assay for transposase-accessible chromatin using high-throughput sequencing (ATAC-seq), which lever ages a robust, transposase enzyme-based method to profile open, accessible chromatin, and requires substantially less input (50,000 nuclei or fewer)^[50]18 and using the remainder for transcriptional analysis by RNA sequencing (RNA-seq), we identified common accessible chromatin regions present in six different mouse organs, as well as the DNA-binding motifs within these regions, to define a ‘core’ endothelial transcriptional code involving the ETS (E26 transformation-specific or erythroblast transformation specific) and SOX (SRY-related HMG box) family of transcription factors. We then identified organ-specific, accessible endothelial gene regulatory regions in embryonic and postnatal development, as well as in the adult mouse. Analysis of these regulatory regions identified transcription factor DNA-binding motifs enriched within distinct organs, while gene expression analysis identified the specific transcription factor family member(s) probably binding these uniquely accessible DNA regulatory elements. Further analysis of the regulatory landscape of the CNS endothelium using single-cell RNA-seq coupled with computational prediction algorithms determined that many of the same gene regulatory networks identified in our ATAC-seq data may govern specific stages of angiogenesis and blood–brain barrier maturation. Finally, we show that transcription factors identified in liver and brain endothelium by our approach are sufficient to molecularly and functionally reprogram naive ECs in vitro. Thus, we present a compendium of shared, and unique, transcriptome and epigenetic data across multiple organs, throughout development, for identification of the transcriptional regulatory networks governing organ-specific endothelial gene expression. Results Using INTACT and ATAC-seq to map accessible chromatin in ECs To analyze organ-specific differences in endothelial chromatin accessibility and gene expression, we combined an endothelial-specific CreERT2 driver line (Cdh5-PAC-CreER)^[51]16 with a Cre-dependent myc epitope-tagged, nuclear envelope GFP reporter mouse (Rosa26^CAG-lox-stop-lox-Sun1-2xsfGFP-6xMyc, or R26^Sun1-sfGFP)15 to enable affinity pulldown of superfolder GFP (sfGFP)-tagged nuclei. Both input (total) and endothelial (GFP^+) nuclei were processed for ATAC-seq^[52]18 and RNA-seq to identify organ-specific accessible chromatin and transcriptional signatures ([53]Fig. 1a). ECs from the embryonic day 12.5 (E12.5) trunk, brain and heart, as well as the postnatal day 6 (P6) and adult mouse brain, retina, heart, lung, liver and kidneys, were analyzed ([54]Supplementary Table 1). Fig. 1 |. Isolation and characterization of tissue-specific endothelial signatures throughout development. Fig. 1 | [55]Open in a new tab a, Workflow for genetic affinity tag labeling of ECs using Cdh5-PAC-CreERT2 and R26^Sun1-sfGFP mice for INTACT. Far left: representative schematic of a blood vessel with GFP-tagged nuclei. Nuclear isolation was followed by RNA-seq profiling of nuclear transcripts and ATAC-seq mapping of accessible chromatin and aligning reads to the mouse genome (far right). b, Far left: various tissues and timepoints used to map endothelial cell diversity in the developing (E12.5), postnatal (P6) and adult (2 months of age) mouse. Far right: representative genome browser tracks from ATAC-seq highlight accessible chromatin regions unique to organ-specific genes such as Map2 in neurons, Rnf207 in cardiomyocytes, Gckr in hepatocytes, Ager in alveolar cells of the lung and Magi-2 in proximal tubule cells of the kidney, and ECs including Cdh5, Pecam1 and Erg. c, Volcano plots show differentially expressed genes (DESEQ2 Wald test two-tailed P values in which differential genes were counted if log[2](FC) > 0.5 and alpha < 0.1.) between input nuclei (blue) and endothelium (red). Note that all developmental timepoints (E12.5–adult) are combined and treated as a single sample for these analyses. To confirm our pipeline, we analyzed chromatin accessibility at loci whose transcripts are enriched in the major (non-endothelial) cellular constituents of each organ sampled. Accessibility at the Map2 ([56]refs. 19,[57]20) locus was enriched in brain input compared with ECs, while Rnf207 ([58]refs. 21,[59]22) was elevated in the heart input, Gckr^[60]23,[61]24 in the liver input, Ager^[62]25 in the lung input and Magi-2 ([63]refs. 26,[64]27) in the kidney input. Conversely, chromatin accessibility around the pan-endothelial markers Cdh5 (encoding vascular endothelial cadherin (VE-cadherin))^[65]28, Pecam1 (encoding CD31)^[66]29 and Erg (encoding ERG)^[67]9 was increased in EC nuclei ([68]Fig. 1b). RNA-seq analysis confirmed organ purity, Sorcs3 ([69]ref. 30) and Map2 were enriched in the brain, while Rnf207 was expressed in the heart^[70]31, Gckr in the liver^[71]32, Ager in the lung^[72]33 and Magi-2 in the kidney^[73]34. These same transcripts were depleted in EC nuclei of their respective organ. Conversely, Mfsd2a was elevated in brain ECs^[74]35, Wt1 in cardiac ECs^[75]36,[76]37, Oit3 in liver ECs^[77]38, Car4 in lung ECs^[78]10 and Gata5 in kidney ECs^[79]39. Pecam1 and Cdh5 were enriched in the endothelial nuclei of all organs ([80]Fig. 1c). We detected little enrichment of transcripts, or chromatin accessibility, for lymphatic EC markers ([81]Extended Data Fig. 1a,[82]b). ECs share core epigenetic features across organs and time We identified 2,646 accessible regions enriched in the endothelium of all organs ([83]Fig. 2a and [84]Supplementary Table 2). We then used Genomic Regions Enrichment of Annotations Tool (GREAT)^[85]40 to link these open chromatin regions to genes to perform Gene Ontology (GO) analysis. Vascular development, blood vessel morphogenesis and angiogenesis were among the top GO terms common in all ECs ([86]Fig. 2b,[87]c and [88]Supplementary Table 3). Hypergeometric Optimization of Motif Enrichment (HOMER)^[89]41 identification of transcription factor DNA-binding motifs enriched in these accessible chromatin regions determined that binding sites for the ETS family of transcription factors (TFs) (ERG, FLI1 and so on), which are crucial for endothelial development^[90]42–[91]45, were the most enriched TF motifs in these common accessible regions ([92]Fig. 2d and [93]Supplementary Tables 4–[94]6). Notably, validated ETS-dependent enhancers in Dll4 (encoding delta-like 4)^[95]44 and Eng (endoglin)^[96]46 were both identified by our analyses ([97]Fig. 2e,[98]f). Fig. 2 |. ECs from diverse organs share a core epigenetic signature. Fig. 2 | [99]Open in a new tab a, Venn diagram showing the overlap of open chromatin regions (2,646 peaks) between the mouse heart, kidney, liver, lung, retina and brain endothelium. b, Heatmap of shared peaks across the endothelial and input datasets. c, GREAT analysis of common peaks showing GO terms related to cardiovascular development and angiogenesis, among others. d, Top 20 transcription factor DNA-binding motifs in shared peaks along with their P value as determined using HOMER (one-tailed Fisher’s exact test). e, Top: position weight matrix (PWM) for transcription factor DNA binding sites of two ETS family members, with the inset box showing the frequency of motif occurrence in all endothelial cells (shared EC) as the distance from the center of the peak within accessible DNA regions as determined by ATAC-seq. f, Representative genome browser tracks from ATAC-seq data highlight accessible chromatin regions in the adult endothelium for the loci surrounding the genes Dll4 and Eng, with putative DNA binding sites for ERG and FLI highlighted in vertical yellow bars and red rectangles. Motifs for the SOX family were the second most abundant known DNA-binding sites enriched in the accessible chromatin of ECs ([100]Fig. 2d). SOXF subfamily members (Sox7, Sox17 and Sox18) have partial redundancy in regulating angiogenesis and vascular maintenance^[101]47–[102]49, and Sox17 regulates arterial differentiation in mice^[103]50. The SOXB1 subfamily member Sox2 has been implicated in endothelial differentiation in vitro^[104]51 and in cerebral arteriovenous malformation^[105]52. Motifs for the Forkhead Box (FOX) family member FOXO1, which regulates angiogenesis and EC metabolism^[106]53–[107]56, were enriched across all organs. Notably, transcripts for multiple ETS (ETS1, ERG, FLI), SOX (SOX17, SOX6, SOX7, SOX4) and FOX (FOXO1) family of TFs were detected in the ECs of all organs ([108]Extended Data Fig. 2a). Identification of organotypic EC transcriptional regulators Merging the E12.5, P6.5 and adult timepoints from each organ into a single dataset, we identified 45,075 EC-enriched, differentially accessible chromatin peaks ([109]Fig. 3a and [110]Supplementary Table 7). These regions were linked to nearby genes^[111]40, and these genes were filtered for enriched expression in the endothelium using our nuclear RNA-seq data (q value < 0.1 and log[2](fold change (FC)) > 0.5). These genes were then used to identify enriched biological processes in each organ ([112]Fig. 3b and [113]Supplementary Tables 8 and [114]9). Brain-enriched accessible chromatin regions were associated with the WNT (wingless-type mouse mammary tumor virus integration site) signaling pathway, while the vasculature of the heart and kidney showed enrichment of semaphorin–plexin signaling, and the heart showed enrichment for Notch signaling. Fig. 3 |. Profiling accessible chromatin and expressed transcripts identifies organ-specific endothelial signatures. Fig. 3 | [115]Open in a new tab a, Heatmap showing differentially accessible regions of open chromatin in the mouse brain, retina, heart, liver, lung and kidney endothelium (45,075 peaks) identified by ATAC-seq. b, Top biological processes from GREAT analysis across differentially accessible peaks in each organ. Only regions annotated to endothelial-enriched genes (determined by RNA-seq) were used in the analysis. c, Top transcription factor motifs in regions of open chromatin in the brain and heart (one-tailed Fisher’s exact test). d, Enriched motifs found by HOMER in each organ. PWM shown over frequency of motif as the distance from the peak center. e, Representative genome browser tracks from ATAC-seq highlighting accessible organ-specific chromatin regions in endothelial-enriched transcripts. Vertical gray bars indicate unique, organ-specific peaks shared across organs with an ERG or FLI DNA binding motif. Comparison of differentially enriched TF motifs across all organs, and filtering by TF expression, revealed that DNA-binding sites for WNT signaling-related transcription factors (ZIC3, LEF1, TCF3 and so on) were enriched in the CNS endothelium ([116]Fig. 3c and [117]Extended Data Fig. 2b,[118]c). Canonical WNT signaling is a key regulator of blood–brain barrier (BBB) development^[119]57–[120]60. Motifs for FOX transcription factors were also overrepresented in the CNS endothelium. While Foxf2 has been implicated in regulating BBB-related genes^[121]60,[122]61, and Foxl2 transcripts are enriched in brain ECs^[123]60 ([124]Extended Data Fig. 2b,[125]c), a requirement for FOXP1, FOXK1, FOXF1 and FOXA1 in BBB function has not been reported. Motifs for the zinc family of transcription factors (GATA1, GATA2, GATA4 and GATA6) were enriched in the heart and liver ([126]Fig. 3c and [127]Extended Data Fig. 2b). GATA1 may regulate EC function in the heart and liver^[128]62, while GATA2 regulates vascular integrity in mice^[129]63. GATA4 and GATA6 expression was enriched in liver and heart ECs and are both required for cardiovascular development^[130]64–[131]66 and endothelial maturation in the liver^[132]67,[133]68. Motifs for nuclear factor of activated T cell (NFAT) TFs were enriched in the cardiac endothelium ([134]Fig. 3c and [135]Extended Data Fig. 2b,[136]c). Nfatc1 is a canonical endocardial marker and is required for cardiac morphogenesis^[137]69,[138]70 and coronary vessel angiogenesis^[139]71, while Nfatc3/c4 null embryos die at E11.5 with excessive vascular growth^[140]72. Motifs for helix-turn-helix (HTH) and winged helix Regulatory Factor binding to the X-box (RFX) are also enriched in the heart ([141]Fig. 3a–[142]c and [143]Extended Data Fig. 2b)^[144]73. Transcripts for Rfx1 and Rfx2 were evident in the endothelium of all organs, and HTH-X-box has an established role in the heart^[145]74. SOX2, SOX3 and SOX4 motifs were present in the ECs of all organs and enriched in the heart ([146]Fig. 3c and [147]Extended Data Fig. 2b). While a role for Sox2 and Sox3 in the cardiac endothelium or endocardium has not been shown, Sox4 is required for outflow tract morphogenesis^[148]75. Motifs for the large musculoaponeurotic fibrosarcoma (MAF) basic leucine zipper (bZip transcription factors) MAFA and MAFB were enriched in the liver endothelium, although MAFB was more highly expressed ([149]Extended Data Fig. 2b,[150]c). MAF TFs interact with ETS1 and SOX TFs^[151]76, and MAFB is involved in angiogenesis^[152]77. Notably, c-MAF was not present in our motif analysis^[153]78,[154]79. These DNA-binding motifs were enriched in the center of regions of open chromatin ([155]Fig. 3d and [156]Supplementary Tables 10–[157]14), suggesting that these factors may drive chromatin accessibility. Several of these unique, tissue-specific accessible regions occurred within, or near, loci of transcripts that are elevated in these organs, and also featured DNA-binding sites for TFs identified in our analyses ([158]Fig. 3e). For example, Slc7a1(encoding solute carrier family 7 member 1), which is enriched in the adult cerebrovascular endothelium^[159]80,[160]81, contains a region downstream from the transcriptional start site (TSS) that is uniquely accessible in CNS ECs and also contains a consensus LEF1 motif ([161]Fig. 3e). This same pattern is found in endocardial genes (Cytl1(encoding cytokine-like 1))^[162]82 (NFAT motifs), liver genes (Dpp4(dipeptidyl-peptidase 4))^[163]83 (GATA4 motifs), lung genes (Ace(angiotensin-converting enzyme)) (FOXO3 motif) and kidney genes (Lgr5(leucine-rich repeat-containing G-protein coupled receptor 5))^[164]84 (HOXC9 motif). Chromatin accessibility change during EC maturation We identified 22,182 peaks specific to the CNS endothelium that were differentially regulated from E12.5 to adult ECs. After annotating peaks to nearby genes, we filtered these data to select for genes whose transcripts were enriched in the endothelium in our RNA-seq data (q value < 0.1 and log[2](FC) > 0.5) for GO analysis (false discovery rate (FDR) < 0.05) ([165]Fig. 4a,[166]b and [167]Supplementary Tables 15–[168]18). We then rank-ordered enriched DNA-binding motifs at each timepoint by the actual expression of each TF (1 = highest, 20 = lowest) ([169]Fig. 4c). At E12.5, FOXF1 and SOX17 were the most highly expressed TFs of those motifs that were enriched in the cerebrovasculature. Members of the ETS family (ETS1, ERG and FLI1) rounded out the remaining top 5 TFs. ETS, FOX and SOX family members made up the rest of the top 20 TFs, as well as TEAD1 and JUN. At P6, ETS1 and FLI1 were decreased, while EHF, ELF3, ELF5 and ERG, along with FOXL2, comprised the top 5. In the adult CNS endothelium, FOXL2 and FOXF1 were followed by ETS1, FLI1 and LEF1. LEF1 and TCF3, regulators of canonical WNT signaling, as well as PPARA, FOXP1, FOXO1, FOXM1, KLF1, KLF5 and NR2F1, were among the other adult-enriched TFs. The same rank ordering of TF motif abundance filtered by mRNA expression was performed for the remaining organs ([170]Extended Data Figs. 3–[171]6 and [172]Supplementary Tables 19–[173]35). Examination of established organ-specific transcriptional regulators, such as Nfatc1 (heart), Gata4 (liver) and Tcf7l2 (brain), revealed differentially accessible potential gene regulatory elements within these loci ([174]Extended Data Fig. 7). Fig. 4 |. Chromatin accessibility changes across time in the brain endothelium. Fig. 4 | [175]Open in a new tab a, Differential chromatin accessibility determined by ATAC-seq within the brain and retinal endothelium (6,540 peaks) of E12.5, P6 and adult mice. b, Top biological processes from GREAT analysis across differentially accessible peaks at each timepoint. c, Top 20 transcription factors ranked by expression for each age. The log[2] expression over input is indicated in the y axis. The motif enrichment P value is shown according to the size of the bubble (one-tailed Fisher’s exact test). d, Normalized gene expression in either E12.5, P6 or adult brain and retina endothelium (top) and genomic tracks for endothelial and input brain samples for genes upregulated in E12.5 (Adm), P6 (Tnc) or adult (Slc9a2). Unique peaks to those timepoints are indicated by the transparent vertical gray bar, and DNA binding sites of the top 20 transcription factor motifs that are present in such peaks are indicated below. avg, average. Examination of brain-specific regions of open chromatin associated with genes that were differentially expressed in E12.5, P6 or adult CNS endothelium showed that Adm (encoding adrenomedullin), enriched in brain EC tip cells^[176]85, contains an accessible chromatin region in E12.5 at the zenith of Adm expression, but this peak is decreased at P6 and not detected in the adult. Similarly, expression of Tnc(encoding tenascin-c), whose protein product is involved in cell adhesion^[177]86, peaks at P6 and features two regions of open chromatin at this same stage that are not detected in the adult endothelium. Finally, Slc9a2, which encodes a Na^+/H^+ exchanger present in brain endothelium^[178]87, contains three regions of open chromatin upstream of its promoter that are enriched in the adult CNS endothelium. Each of these accessible chromatin regions contain putative binding sites for the relevant TFs identified (for example, ETS1, FOXF1 and FOXO1 at E12.5 in Adm; TCF3, FOXF1 and LEF1 at the adult in Slc9a2) ([179]Fig. 4d). EC transcription factor motifs are evolutionarily conserved To determine whether the transcriptional networks identified in the mouse brain are conserved in humans, we turned to an in vitro model of the human brain endothelium: hCMEC/D3 cells^[180]88–[181]91. A total of 93,197 open chromatin regions in human brain microvascular ECs and 13,241 peaks enriched in P6 and adult mouse brain ECs were used to identify enriched motifs. DNA-binding sites for ERG, FOXF1, FOXO3, FOXL2, FOXP1, TCF3, LEF1, FOXO1, TCF7 and ZIC3 were enriched in accessible chromatin of the adult human and mouse cerebral endothelium ([182]Fig. 5a,[183]b). These conserved regions of open chromatin were associated with processes such as vascular remodeling, cell junction assembly and negative regulation of the cell cycle (via BMP/SMAD and vascular endothelial growth factor (VEGF) inhibition) ([184]Fig. 5c and [185]Supplementary Tables 35 and [186]36). Examples of evolutionarily conserved, non-coding regions containing these motifs are found in Slc31a1 and proximal to Mfsd2a ([187]Fig. 5d,[188]e). Fig. 5 |. Evolutionary conservation of transcription factor motifs in regions of open chromatin between humans and adult mice. Fig. 5 | [189]Open in a new tab a, Diagram representing the total number of open chromatin regions in hCMEC/D3 cells and the adult endothelial brain vasculature. b, Selected known transcription factor DNA binding motifs in conserved peaks along with their P value after analysis by HOMER (one-tailed Fisher’s exact test). c, GO term analysis of genes with conserved nearby accessible chromatin regions. d,e, Representative genome browser tracks of Slc31a (d) and Mfsd2a (e) highlighting (in yellow) conserved accessible chromatin regions in humans (top) and mice (bottom) as defined by ATAC-seq and Omni-ATAC-seq. Transcription factor motifs present in the highlighted peak are shown above. Conservation at the nucleotide level within each highlighted peak is shown below each locus. Cons, conservation; verts, vertebrate. Single-cell resolution of BBB maturation We next used single-cell RNA-seq (scRNA-seq) and in silico prediction of endothelial gene regulatory networks as an orthogonal approach to validate our ATAC-seq and RNA-seq results. CD31^+ ECs from whole brains (E9.5, E12.5 and E16.5) or the cortex (P8 and adult) were isolated by Magnetic Activated Cell Sorting (MACS) and processed for scRNA-seq. After filtering ([190]Methods), 22,990 cells were subjected to dimensionality reduction analysis for visualization. Assigning cell identities using established markers revealed the presence of microglia, mural cells, astrocytes or radial glia, and endothelium ([191]Supplementary Table 37 and [192]Extended Data Fig. 8a–[193]d). We used NicheNET^[194]92 to predict cell–cell-based interactions between components of the neurovascular unit (NVU) from these scRNA-seq data^[195]93. Designating ECs as signal-receiving cells (receptors and downstream effectors) and all other cell types of the NVU (neurons, mural cells, glia, astrocytes) as signal-sending cells, we examined signaling networks that were upregulated in the adult brain compared with the E9.5 brain ([196]Extended Data Fig. 8e). NicheNET predicted that pericyte-derived cadherin-2 (N-cadherin) may activate Ctnnb1 (encoding catenin beta-1/β-catenin) and Ctnnd1(catenin delta 1; P120)^[197]94. It also predicted that pericyte-presented CDH2 (cadherin-2, N-cadherin) activated WNT-related genes in ECs, including Cdkn1a^[198]95, Ccnd1 ([199]refs. 96,[200]97) and Ptma^[201]98. While most activation of Ctnnb1 in the developing CNS occurs via neural-precursor-cell-mediated secretion of WNT ligands, N-cadherin^+ pericytes can adhere to ECs expressing a different cadherin (that is, VE-cadherin), and binding of pericyte N-cadherin to EC-presented cadherins could impact the cytoplasmic pools of P120 and β-catenin within ECs^[202]99,[203]100. Furthermore, pericyte-specific Cdh2 deletion increased permeability and reduced VE-cadherin-positive adherens junctions in neighboring ECs^[204]101,[205]102, and pericyte CDH2 induces expression of WNT target genes in ECs^[206]103,[207]104. NicheNET also predicted that smooth muscle cell-presented Itgb1 (encoding integrin beta 1) may activate Rad51, which is involved in vascular maintenance^[208]105,[209]106. We next used CCInx (version 0.5)^[210]107 to identify differential ligand–receptor signaling between the endothelium and pericytes in the young versus adult mouse brain ([211]Extended Data Fig. 8f). CCInx predicted that pericyte vitronectin interacted with endothelial integrin β[1] in adult mice. Recent studies showed that pericyte-derived vitronectin-mediated activation of integrin α[5]β[1] in ECs is essential for BBB function^[212]108, and integrin β[1] is required for BBB function^[213]109–[214]111. CCInx also suggested that the chemokines that regulate leukocyte migration and maintain homeostasis, such as Cxcl12, are enriched in adult brain ECs, while its receptor, Ackr3/Cxcr7, is enriched in pericytes^[215]112,[216]113. Similarly, adult brain ECs express Pdgfb, while its cognate receptor, Pdgfrb, was enriched in adult pericytes^[217]114,[218]115. Conversely, Mif(encoding macrophage migration inhibitory factor), which attenuates pericyte contractility^[219]116, was reduced in adult brain ECs, and expression of multiple potential MIF receptors, including Tfrc/Cd71 (transferrin receptor 1) and Itgα4(integrin α4), was reduced in pericytes. Thus, cellular communication within the NVU can be inferred from these results, and these EC-enriched datasets may be integrated with other existing scRNA-seq datasets for further analyses. In total, 22,990 brain ECs (Pecam1^+, Cdh5^+) from each timepoint (E9.5 = 6,039; E12.5 = 6,822; E16.5 = 3,358; P8 = 4,048; adult = 2,723) were extracted and combined for analysis. Dimensionality reduction and visualization^[220]117,[221]118 revealed that while EC-specific transcripts (that is, Pecam1, Cdh5) were present across all timepoints, canonical lymphatic endothelial cell markers (that is, Lyve1, Prox1) were lowly expressed, suggesting limited presence of non-blood ECs ([222]Fig. 6a and [223]Extended Data Fig. 9a,[224]b). Fig. 6 |. Maturation of the BBB at single-cell resolution. Fig. 6 | [225]Open in a new tab a, UMAP of all only subset ECs (from [226]Extended Data Fig. 8) colored coded by timepoint and clustered using Monocle 3 to show maturation across time. b, A feature plot showing the expression score of a 66-marker panel BBB gene signature shows later acquisition of maturing BBB markers. c,d, Feature plots for the BBB maturation markers Mfsd2a (c) and Abcb1a (d). e, Confocal micrographs of sections processed for smFISH for Mfsd2a transcripts (magenta), along with antibody staining for the EC marker Cldn5 (yellow) and DAPI staining for nuclei (turquoise) from E12.5, P7 and adult mouse brains. f, Confocal micrographs of sections processed for smFISH for Abcb1a transcripts (magenta), along with antibody staining for the EC marker Cldn5 (yellow) and DAPI staining for nuclei (turquoise) from E12.5, P7 and adult mouse brains. In e and f, the area in the dashed white box is magnified in the image below. g, SCENIC binary activity heatmap representing active regulons in brain ECs across all timepoints. Vertical columns represent individual sequenced cells, while each horizontal row represents an individual regulon. Highlighted regulons are shown in the right. Differentially active regulon in the cerebral endothelium at E9.5 and E12.5 (green shading) compared with E16.5, P8 and adult (blue shading) and P8 and adult (orange shading). Highlighted for each group are representative regulons with their binding motif visualization logo and the number of genes that comprise the regulon in parentheses. h, Confocal micrographs of sections processed for smFISH for Cldn5 (top), Plagl1 (middle) and Dbp (bottom) transcripts (magenta), along with antibody staining for the EC marker Cldn5 (yellow) and DAPI staining for nuclei (turquoise) from E12.5, P7 and adult mouse brains. In the images, the area in the dashed white box is magnified in the image below. Differential gene expression signatures were evident between the timepoints ([227]Extended Data Fig. 9c). Interestingly, Mif, which encodes an inflammatory cytokine implicated in angiogenesis^[228]119, was robustly expressed in E9.5 brain ECs, but downregulated at later stages. Marcksl1, a gene involved in blood vessel morphogenesis^[229]120, was an upregulated gene in the E12.5 brain endothelium. Expression of Pvlap, a pan-endothelial marker that is lost in the mature BBB endothelium^[230]121,[231]122, decreased after E12.5. Conversely, the BBB markers Mfsd2a^[232]35,[233]123 and Slco1a4 ([234]refs. 124,[235]125) were highly expressed at E16.5 or later, coinciding with BBB maturation^[236]35,[237]60. Overall, BBB-related genes^[238]126 ([239]Supplementary Table 38) were enriched at later timepoints ([240]Fig. 6b–[241]d and [242]Extended Data Fig. 9d). Immunofluorescence for the EC marker CLDN5, combined with single-molecule fluorescence in situ hybridization (smFISH) for Mfsd2a and Abcb1a ([243]Fig. 6e,[244]f), confirmed the trends shown in our scRNA-seq analyses. Gene regulatory networks regulating brain EC maturation We next correlated expression of TFs and their presumptive transcriptional targets to identify differentially active gene regulatory networks (GRNs) during cerebrovascular maturation. Using Single-Cell Regulatory Network Inference and Clustering (SCENIC)^[245]127, we identified three distinct clusters of ECs ([246]Fig. 6g). The first cluster was enriched primarily in E9.5 and E12.5 brain ECs and included GRNs (or regulons) driven by SOX11 (106 genes), PLAGL1 (79 genes) and E2F1 (314 genes). SOX11 regulates vascular development and is active during pathological angiogenesis^[247]128,[248]129, while PLAGL1 may impact early developmental angiogenesis^[249]130, and E2F1 modulates VEGF expression in vitro^[250]131. The second cluster of regulons were active primarily in the E16.5, P8 and adult CNS endothelium, and included JUN (33 genes), FOXF1 (110 genes) and FOXQ1 (236 genes). Jun has been implicated in sprouting angiogenesis^[251]132–[252]134, while Foxf1 is required for EC barrier function in the lung^[253]135. Foxq1 is enriched in the developing brain endothelium, but is not sufficient to induce expression of BBB markers^[254]60. The third major cluster, enriched in the P8 and adult CNS endothelium, contained LEF1 (43 genes), DBP (96 genes), FOXP1 (532 genes) and NR3C1 (369 genes) GRNs. Lef1, which encodes an obligate binding partner of β-catenin in the nucleus, acts downstream of canonical WNT signaling to govern BBB function^[255]136–[256]142. Nr3c1(nuclear receptor subfamily 3 group C member 1), which encodes a glucocorticoid receptor, regulates the WNT–β-catenin pathway^[257]143. Notably, Dbp(encoding albumin D-binding protein), a proline amino-acid-rich domain bZip (PAR bZip) transcription factor involved in circadian rhythm control in the BBB^[258]144,[259]145, was increased in adult CNS ECs. smFISH confirmed that expression of BBB components, such as Cldn5 (encoding claudin-5), increases over time. This analysis also showed that transcriptional regulators identified by SCENIC, such as Plagl1, were enriched in the vasculature of the embryonic brain, whereas adult-enriched transcriptional regulators such as Dbp were robustly detected in the vessels of the P7 and adult mouse brain ([260]Fig. 6g,[261]h). Notably, many of the same GRNs identified by SCENIC (JUN, FOXF1, LEF1 and so on) were initially identified by our ATAC-seq approach ([262]Fig. 4c). Developmental trajectories of the mouse cerebrovasculature Our scRNA-seq data contained a mixture of mitotic (Top2a, Mki67, Aurka), arterial (Fabp4, Fbln5, Vegfc, Gkn3, Bmx), angiogenic (Meox1, Angpt2), capillary (Hist1h1b, Car14, Gpatch1) and venous (Fmo2, Il1r1, Nr2f2, Ackr1) ECs ([263]Fig. 7a–[264]c and [265]Supplementary Table 39). GO analysis revealed that Cap-1 and mitotic cells were enriched for processes related to mRNA splicing, cell cycle and mitosis. Mature clusters, such as Cap-2 and Art-2, were enriched for extracellular matrix organization, cell adhesion and negative regulation of cell proliferation ([266]Extended Data Fig. 10). Fig. 7 |. Temporal transcriptomic profiling of the mouse cerebrovasculature. Fig. 7 | [267]Open in a new tab a, UMAP of ECs (from [268]Fig. 6) colored coded by timepoint and clustered using Monocle 3 to show maturation across time. b, The same trajectory plot color coded by EC subtype. c, Feature plots for pan-EC (Cdh5, Pecam1), mitotic (Top2a), Cap-1 (Car2), Cap-2 (Fmo2), venous (Il1r1), angiogenic (Meox1), Art-1 (Igf1) and Art-2 (Hey1) markers. d, Pseudotime trajectory (solid black line; the dashed arrow indicates the direction of pseudotime progression) of arterial maturation projected onto the UMAP. Pseudotime scores are plotted for each single cell across the UMAP embedding in the figure. Gray cells are excluded and have no value. e, Arterial gene expression trends across pseudotime. f, Confocal micrographs of sections processed for smFISH for Stmn2 transcripts (magenta), along with antibody staining for the EC marker Cldn5 (yellow) and the smooth muscle marker SMA (white) to highlight arteries, and DAPI staining for nuclei (turquoise), from E12.5, P7 and adult mouse brains. The areas in the white boxes are magnified in the images below. g, Psuedotime trajectory of venous maturation. h, Venous gene expression trends across pseudotime. i, Confocal micrographs of sections processed for smFISH for Foxq1 transcripts (magenta), along with antibody staining for the EC marker Cldn5 (yellow) and the smooth muscle marker SMA (white) to distinguish arteries from veins, and DAPI staining for nuclei (turquoise), from E12.5, P7 and adult mouse brains. The areas in the white boxes are magnified in the images below. j, Heatmap showing the expression of arterial and venous gene modules across various populations. k, Feature plots of module expression across endothelial maturation. High value scores (red) indicate that the genes within the module are expressed higher than a set of randomized control genes in a single cell. Low module scores (blue) indicate that the genes in that module are expressed lower in aggregate compared with the background. A score of zero represents a module with equal expression to background. Given the complexity of endothelial cell ontogeny, we separated arterial-like ECs for pseudotime trajectory analysis^[269]146. This analysis revealed a progression from a naive, mitotic identity to an intermediate arterial identity (Art-1), and ended at a mature arterial state (Art-2) ([270]Fig. 7d). Markers of arterial maturation reflected this progression. For example, the expression of Plvap, which regulates selective permeability of endothelial fenestrations^[271]147, is elevated during development but decreased during BBB formation ([272]Fig. 7e). Conversely, Stmn2 (stathmin-2), which encodes a microtubule-associated protein, was enriched in mature arterial endothelial ECs. Whether Stmn2 is a true marker of arterial maturation was unknown^[273]148,[274]149. smFISH to visualize Stmn2 transcripts, combined with counterstaining for CLDN5 to label ECs, or smooth muscle actin staining to differentiate arteries from veins, confirmed our scRNA-seq results and validated Stmn2 as a marker of maturing arterial ECs ([275]Fig. 7f). The venous endothelium initiated from a mitotic origin in E9.5 and E12.5 ECs, and then transitioned through an early capillary (Cap-1) to an angiogenic (Angio) phase, and terminated in either a mature capillary (Cap-2) or venous (Venous) identity ([276]Fig. 7g). Within this trajectory, the Apln expression peaked and then sharply decreased in mature venous ECs ([277]Fig. 7h), consistent with previous reports^[278]150. smFISH showed that the expression of Foxq1, a transcription factor identified by SCENIC, increased in maturing venous ECs ([279]Fig. 7i). Overall, these vignettes show that the trajectories derived from our scRNA-seq studies accurately inform on changes in endothelial cerebrovascular maturation. We next quantified gene expression changes as a function of pseudotime (q value < 1 × 10^−100). In total, 2,197 genes were differentially expressed across the arterial trajectory, while 2,866 transcripts changed during venous maturation. To avoid the pitfalls of relying only on clustering-based approaches^[280]151, we merged genes with similar expression patterns into co-expression modules ([281]Fig. 7j). The arterial time course produced 16 modules, while the venous trajectory consisted of 20 modules ([282]Supplementary Tables 40 and [283]41). While some modules overlapped between the trajectories (for example, arterial module 3, venous module 19), many were specific to one fate (for example, arterial module 7, venous module 5). Some modules were cluster specific, such as venous module 17, which marked angiogenic ECs. We next performed pathway enrichment analysis to predict the biological function(s) of these modules. Many of the enriched pathways in the modules were consistent with the pseudotime analyses, as arterial modules 1, 9, 12 and 14 showed enrichment for the cell cycle ([284]Supplementary Fig. 1a and [285]Supplementary Tables 42 and [286]43), and these same modules were enriched in the early mitotic cluster ([287]Fig. 7b,[288]k). Similarly, modules 1, 4 and 19, which were present in the mature Cap-2 and venous clusters, were enriched for blood vessel morphogenesis and positive regulation of cell migration ([289]Supplementary Fig. 1b,[290]c and [291]Supplementary Tables 44 and [292]45). Given the phenotypic plasticity of the endothelium^[293]152, we hypothesized that the expression of TFs enriched in mature CNS ECs should reprogram the endothelium to adopt a more BBB-like identity, while transcriptional regulators enriched in the leaky, sinusoidal endothelium of the liver should have the opposite effect. To determine the optimal combination of BBB TFs, we focused our studies on those regulators identified by our bulk RNA-seq and ATAC-seq studies that were also identified in our scRNA-seq analyses, specifically FOXF1, FOXF2, LEF1, ZIC3, TCF7, TCF7L1 and TCF7L2. Expression of these individual factors increased electrical resistance across an endothelial monolayer (transendothelial electrical resistance (TEER)) compared with control cells. Conversely, a dominant negative TCF4 (dnTCF) had the opposite effect. At 48 h, TCF7 and FOXF1 yielded the highest resistance of the individual TFs. Given previous data implicating FOXF2 ([294]ref. 60), we assayed their combined impact and found that they increased resistance compared with either factor alone (P < 0.05). These data were confirmed using Transwell permeability assays ([295]Supplementary Fig. 2a,[296]b). We then tested various TF combinations to identify the optimal BBB-inducing TF cocktail. LEF1/ZIC3/TCF7/FOXF1/FOXF2, followed by LEF1/ZIC3/TCF7L2/FOX1/FOX2, enhanced barrier function in the TEER and Transwell assays over all other combinations, and provided the most robust induction of BBB-associated channels and transporters (ABCB1A, SLCO1A4, SLCO1C1, MFSD2A) and junctional markers (CLDN5, CDH5) ([297]Supplementary Figs. 2c,[298]d and [299]3a–[300]d). We compared the impact of this optimal BBB TF cocktail to regulators identified in our bulk RNA-seq and ATAC-seq mouse liver data, specifically MAFB, GATA2 and GATA4. While a role for MAFB in barrier function and an interaction with GATA2 and GATA4 in ECs have not been shown, our bulk RNA-seq and ATAC-seq studies revealed MAFA and MAFB DNA-binding motif enrichment (but only elevated Mafb transcripts). In addition, GATA2 and GATA4 are required for development^[301]153,[302]154, and GATA4 overexpression represses continuous EC genes (Pecam1, Cav1, Emcn, Acvrl1, Sox18, Vegfa)^[303]67. Transduction of human umbilical vein endothelial cells (HUVECs) with the BBB factors induced robust impedance, yet the liver TFs, as well as dnTCF4, decreased barrier function ([304]Fig. 8a–[305]c). Markers of BBB maturation (ABCB1A, SLCO1A4, SLCO1C1, MFSD2A) and transcripts encoding junctional proteins (that is, CLDN5) were upregulated by the BBB cocktail, and unaffected or downregulated by dnTCF4 and the liver TFs. Notably, only the liver TFs induced the sinusoidal endothelial markers STAB1 and DNase1L3 ([306]Fig. 8d,[307]f). Immunostaining confirmed BBB TF induction of the tight junction markers CLDN5 and VE-cadherin ([308]Fig. 8g). These functional studies highlight the utility of this resource to identify transcriptional regulators of organospecific EC maturation. Fig. 8 |. Organ-specific transcription factors reprogram endothelial cell gene expression and alter barrier functionality. Fig. 8 | [309]Open in a new tab a, Experimental workflow for HUVEC transduction with organ-specific transcription factor cocktails and subsequent functional analyses. b, Quantification of TEER across an endothelial monolayer over time, with T[0] denoting the day of lentiviral-mediated transduction of the indicated transcription factors. c, Quantification of FITC-conjugated 40 kDa dextran leak across an endothelial monolayer in a Transwell permeability assay following lentiviral transduction with the indicated transcription factors. d–f, qRT–PCR analysis of the induction of BBB-related endothelial markers (d), genes encoding junctional proteins (e) and sinusoidal liver endothelial transcripts (f) following transduction of HUVECs with the indicated transcription factors. Expression is relative to non-transduced cells. Non-trans, non-transduced. g, Representative images of HUVECs immunostained for the junctional proteins VE-cadherin (cyan) and Claudin5 (yellow), with nuclei labeled by DAPI (gray) staining, and visualization of the mCherry reporter (magenta) following transduction with the indicated transcription factor(s). There were three biological replicates, each containing three technical replicates, for each experiment. Data shown are mean ± s.d. *P < 0.05; **P < 0.01; ***P < 0.001; ****P < 0.0001 (two-way ANOVA with Tukey’s multiple-comparison test). Panel a created with [310]BioRender.com. Discussion Herein we have established a lexicon of common, accessible regions of open chromatin present within the endothelium of six organs, across development, and identified a core set of enriched transcription factor DNA binding motifs present in most ECs. Furthermore, we have also identified uniquely accessible chromatin regions in ECs of specific organs, along with the possible transcription factors that act on these regions, which may generate the functional heterogeneity evident within these different vascular beds^[311]85. We also show that these regulatory regions, and the transcription factor motifs within them, are conserved in the human CNS endothelium. Finally, using single-cell transcriptomic approaches, we uncovered numerous novel gene regulatory networks governing cerebrovasculature maturation and validated these putative regulators in vitro. As expected given the role of ETS TFs in endothelial biology, the DNA-binding motif for this family of TFs was the most commonly occurring DNA binding site within accessible chromatin of the endothelium (regardless of organ identity)^[312]155. Other common, core motifs were those of the SOX transcription factor family^[313]48,[314]51. Critically, organ-specific signatures also emerged, such as NFAT, AP1, SOX4 and 3, and RFX1 and 2 and FOXD3 in the heart. By analyzing samples from developing and mature organs, we identified transitions in motif usage and TF expression, such as the shift from ETS factors in angiogenic E12.5 brain ECs to FOX family members in the adult CNS endothelium. This catalog of chromatin landscapes and gene expression signatures of ECs from different organs, across development, can be used to generate new hypotheses regarding endothelial specialization, maturation and homeostasis. Notably, the major DNA-binding motifs found in the mouse adult cerebrovasculature were present in a human cell culture BBB model and shared similar GO signatures. Future studies will interrogate the role these conserved brain-specific regulatory elements play in modulating gene expression in vivo. Given the differences between the pre- and post-BBB endothelium^[315]156, we sought to identify the transcriptional networks mediating CNS maturation. Our scRNA-seq studies identified novel GRNs in the early brain, such as SOX11, PLAG1 and E2F1, and confirmed our ATAC-seq results. Whether the unique TFs identified herein, such as FOXF1, regulate BBB maturation remains unknown. We also identified LEF1, NR3C1 and DBP regulons in the adult brain endothelium. These data are consistent with a requirement for Lef1 in BBB maturation^[316]59,[317]142,[318]157,[319]158, reports suggesting Dbp controls circadian rhythms within the CNS^[320]145,[321]159 and studies suggesting NR3C1 plays a role in aneurysms^[322]160,[323]161. Previous data show FOXF2, FOXL2 and ZIC3 are enriched in the E14.5 mouse brain endothelium, and the expression of FOXF2 and ZIC3 induced the BBB ABCB1, TNFRSF19 and SLCO2B1 in HUVECs^[324]60. Our bulk RNA-seq and ATAC-seq data confirm that these TFs, as well as their DNA-binding motifs, are present in maturing cerebrovascular ECs. However, filtering these DNA-binding motifs against actual TF expression showed that FOXF1 (not FOXF2) was enriched in brain ECs, while scRNA-seq identified FOXF1 GRNs in the adult CNS endothelium. Our functional validation in HUVECs also confirmed that FOXF1 increased barrier function and BBB gene expression more robustly than FOXF2. The four mammalian TCF genes (LEF1, TCF7, TCF7L1 (TCF3) and TCF7L2 (TCF4)) recognize nearly identical DNA sequences (5′-C/GATCAAGC/G-3′), overlap in their chromatin occupancy and show partial functionally redundancy^[325]162–[326]165. While the expression of each individual TCF increased barrier function and BBB-related markers, TCF7 showed the most potent interaction when combined with LEF1, ZIC3, FOXF1 and FOXF2. Conversely, the transduction of TFs enriched in liver ECs (MAFB, GATA2 and GATA4) increased barrier permeability and induced the expression of the sinusoidal EC markers. Future investigations will determine whether TFs enriched during organogenesis or EC homeostasis (that is, E12.5, P6 and adult) control organotypic endothelial maturation. Cell trajectory analysis identified unique gene expression modules within EC subpopulations during the transition from a mitotic and proliferative signature in the early cerebral endothelium to mature arterial, capillary and venous fates in the adult brain. GRNs having differential activity in early angiogenic ECs versus homeostatic adult brain ECs may play a role in controlling endothelial maturation. Future studies will determine whether human BBB differentiation features similar GRNs and whether these trajectories and the GRNs enriched within these populations are altered in animal models of neurovascular disease featuring BBB disruption. Limitations of the present study The transcriptional status of many loci is not determined by differential chromatin accessibility, but by differential recruitment of the transcriptional machinery or posttranslation modification of already engaged protein complexes (for example, Pol II pause–release)^[327]43,[328]166–[329]168. Herein we have linked potential cis regulatory elements to nearby target genes, not directly measured physical contacts, and these putative enhancers have not been validated. Bulk nuclear RNA-seq yielded fewer transcripts than bulk whole-cell RNA-seq. Whether this was due to the loss of cytoplasmic RNA or inadequate input material is unclear. Finally, the chromatin accessibility data from cultured ECs do not reflect the transcriptional complexity of the adult human brain. Methods Mice All mouse protocols were approved by the Institutional Animal Care and Use Committee at Baylor College of Medicine and the University of Virginia. For all experiments, noon on the day a vaginal plug was discovered was considered E0.5, the day of birth was considered P0 and all adult mice were 8 weeks of age. Mice were housed with access to food (normal chow diet) and water ad libitum on a 12-h light–12-h dark cycle at 21 °C and 50–60% humidity. Genotyping and mice used Cdh5-BAC-CreERT2 mice (MGI, number 3848982) were from R. Adams (Max Planck Institute for Molecular Biomedicine). Rosa26-Sun1-sfGFP-6xMyc (for example, R26^Sun1GFP; MGI, number 5443817) was purchased from Jax. Genotyping for all alleles was performed by PCR using the gene-specific primers listed in [330]Supplementary Table 46. Mouse endothelial nuclear isolations For embryonic analysis, tamoxifen (0.015 mg per kg (body weight)) was administered to pregnant dams by intraperitoneal injection at E10.5 and embryos were collected at E12.5. For postnatal tissues, tamoxifen (0.015 mg per kg (body weight)) was administered by subcutaneous injection at P1 and P3, and tissues were collected at P7. For adult experiments, tamoxifen (0.015 mg per kg (body weight)) was administered by intraperitoneal injection 7 days before tissue isolation. In all cases, after gross dissection, GFP expression within the vasculature of each tissue of interest (or embryo) was confirmed by direct immunofluorescence for each sample collected. GFP-negative samples were not processed further. Nuclear isolation was performed according to [331]ref. 15. Briefly, fresh tissue was collected on ice in Buffer HB++ composed of 0.25 M sucrose, 25 mM KCl, 5 mM MgCl[2], 20 mM tricine–KOH, pH 7.8 with protease inhibitors (Roche/Sigma, catalog number 11873580001), 1 mM DTT (Sigma D0632), 0.15 mM Spermine (Sigma S1141), 0.5 mM Spermidine (Sigma S2501) and RNAse inhibitors (Promega N2611) and immediately dissected and minced into 1-mm by 1-mm portions with curved scissors. Tissue was transferred along 1 ml of HB++ in a chilled Eppendorf tube in ice and homogenized using a Bio-gen Series PRO200 homogenizer. Short bursts of ~5–8 s were done to prevent overheating. Once no large pieces were observed, the tissue was transferred to large-clearance dounce homogenizer ‘A’ (7 ml; Kontes Glass Company) and 4 ml of HB++ was added. Tissue was homogenized with 20 strokes and transferred to a small-clearance homogenizer ‘B’; 320 μl of 5% IGEPAL CA-630++ in HB++ was added, and the tissue was dounced with the tight pestle 20 more times, slowly to avoid creating bubbles and disrupting cell membranes. The homogenate was then strained using a 40-μm cell strainer into a 50-ml conical tube. Subsequently, 5 ml of working solution of 5 volumes of Optiprep solution (Sigma, D1556) and one volume of diluent (150 mM KCl, 30 mM MgCl[2], 120 mM tricine–KOH, pH 7.8 in water) were added, homogenized by inversion and poured into an empty prechilled 30 ml Corex tube. Once all samples were ready, with the use of a pipette aid, the tip was placed just above the bottom surface of the Corex tube, and the sample was slowly underlying with 7.5 ml of the 30% and then 4 ml of the 40% iodixanol++ solutions (diluted with buffer HB). Nuclei were then isolated by density gradient centrifugation with Optiprep density gradient medium. Nuclei were collected from the 30–40% interface and then pre-cleared with Protein-G Dynabeads (Life Technologies, 10003D). A portion of these nuclei was held back for use as input samples. Next, GFP^+ nuclei were immunoprecipitated with an anti-EGFP antibody (Thermo Fisher, number [332]G10362; ABfinity rabbit monoclonal anti-GFP antibody; 0.2 mg ml^−1) for 40 min at 4 °C with gentle agitation, followed by binding to Protein-G Dynabeads (Invitrogen, 10003D) for 20 min at 4 °C to enrich for endothelial cell nuclei. Isolated nuclei were filtered using a 20-μm Celltrics filter (Sysmex, number 04-004-2326). Specific amounts of tissue and yields of nuclei from each tissue are listed below. For INTACT and downstream ATAC-seq and nuclear RNA-seq, tissue was collected from mice on a mixed Cdh5-PAC-CreER; Rosa26^lsl-Sun1-sfGFP (FBV and C57/BL6) genetic background. In most cases, tissues were pooled from multiple animals and used as a single biological sample. The total number of biological samples for ATAC-seq and RNA-seq for each organ are shown in [333]Supplementary Table 1. For adult hearts, 4 hearts (2 male, 2 female) were used per INTACT experiment, with 80% of the tissue processed resulting in a total of 1.07 × 10^6 isolated nuclei. For adult lungs, 2 lungs (1 male, 1 female) per INTACT experiment were used with 60% of the tissue processed and resulting in a total of 1.1 × 10^6 isolated nuclei. Ten adult retinas (3 female, 2 male) were used per INTACT experiment resulting in 50,000 isolated nuclei. One adult brain was used per INTACT experiment (3 total experiments, 2 EC samples were 1 male and 1 female, while the 1 input sample was from the male), resulting in 1.45 × 10^6 isolated nuclei. One adult liver was used per INTACT experiment, with 50% of the tissue processed resulting in 5 × 10^5 isolated nuclei. Four adult kidneys (2 male, 2 female) were used per INTACT experiment with 60% of the tissue processed and resulting in 8.5 × 10^5 isolated nuclei. Eight P7 hearts (3 male, 5 female) were used per INTACT experiment resulting in 9 × 10^5 isolated nuclei. Four P7 lungs (2 male, 2 female) were used per INTACT experiment resulting in 395,000 isolated nuclei. Sixteen P7 retinas (3 male, 5 female) were used per INTACT experiment resulting in 85,000 isolated nuclei. Six P7 brains (3 male, 3 female) were used per INTACT experiment with 60% of the tissue processed and resulting in 1 × 10^6 isolated nuclei. Eight P7 livers (3 male, 5 female) were used per INTACT experiment with 60% of the tissue processed and 6.7 × 10^5 isolated nuclei. Eight P7 kidneys (3 male, 5 female) were used per INTACT experiment resulting in 1.2 × 10^6 isolated nuclei. Five E12.5 trunks (sex undetermined) were used per INTACT experiment resulting in 1.0 × 10^6 isolated nuclei. Five E12.5 brains (sex undetermined) were used per INTACT experiment resulting in 1.5 × 10^5 isolated nuclei. Five E12.5 hearts (sex undetermined) were used per INTACT experiment resulting in 47,000 isolated nuclei. Each isolation was performed at least twice. ATAC-seq Approximately 50,000 bead-bound EGFP^+ and 50,000 input nuclei from mouse tissues were used as input for ATAC-seq. ATAC-seq libraries for mouse ECs were processed as previously described^[334]169, and libraries were generated using the Nextera DNA Sample Preparation Kit (Illumina, FC-121-1030). The quality of purified DNA libraries was checked using an Agilent High Sensitivity DNA kit (Agilent Technologies). Paired-end, 2 × 75 bp sequencing was performed on an Illumina NextSeq 500 instrument. Reads were mapped to the mm10 version of the mouse genome using Bowtie2 with default paired-end settings^[335]170. Mitochondrial reads, reads with a mapping quality (MAPQ) < 10 and reads that did not align with the reference genome were removed using Samtools (version 1.13)^[336]171. Duplicated reads were then removed with Picard MarkDuplicates^[337]172. Peak calling was carried out with MACS2 (callpeak –nomodel –broad) (v2.2.7.1)^[338]173. Diffbind (version 3.2)^[339]174,[340]175 was used to import peaksets (min.overlap= 0.66) into RStudio Server (version 1.4.1717; [341]https://www.rstudio.com) using R (version 4.1)^[342]176. The dba.blacklist function was used to filter out peaks that overlap with the ENCODE blacklist. Differentially accessible regions between the endothelium and the input nuclei of each organ were extracted using DESeq2 (version 1.32.0)^[343]177 with P < 0.5 and >1 fold change difference. Endothelial-enriched peaks from each organ were compared using the mergepeaks function in Homer (version 4.11)^[344]41. Peaks present in all organs were used for analysis in [345]Fig. 2. Peaks present in single organs were used for analysis in [346]Fig. 3 and [347]Extended Data Figs. 3–[348]6. Motif enrichment analysis was conducted with findMotifsGenome and enrichment graphed as previously described^[349]178. Graphs for individual motif distance from peaks were created using annotatePeaks in Homer and presented in an enrichment plot^[350]178. GO analysis was done using GREAT (version 4.0.4)^[351]40. Nuclear RNA-seq In parallel with our ATAC-seq experiments, all remaining bead-bound EGFP^+ nuclei were processed for RNA extraction using the RNeasy Plus Micro kit (Qiagen). Nuclear RNA-seq libraries were constructed with the Stranded RNA-seq Kit with Ribo Erase (Kapa Biosystems, KK8484) with custom Y-shaped adaptors. Paired-end 2 × 75 bp NSQ 500/550 Hi Output KT v2.5 – 75 CYS (Illumina, 20024906) was performed for RNA-seq libraries on an Illumina NextSeq 500 instrument. Reads were first mapped to the mouse genome (mm10) using Salmon (version 1.5.1)^[352]179. Transcript-level quantification was then imported using txtimport (version 1.20.0)^[353]180 and analyzed using DESeq2 ([354]ref. 177). Differentially expressed genes between the endothelial and input nuclei were defined as those transcripts with an expression log[2](FC) > 0.5 and Benjamini–Hochberg adjusted P value (q value) < 0.1. Volcano plots were created using EnhancedVolcano (version 1.10.0)^[355]181. Magnetic activated cell sorting for mouse scRNA-seq Brain tissue was processed for CD31 MACS with slight variations depending on the timepoint analyzed. For embryonic brains (E9.5, E12.5, E16.5), embryos from timed C57BL/6 mouse matings were collected in ice-cold buffer HBSS++ (HBSS plus FBS, pen–strep and HEPES). Dissected brains were pooled and then placed in 250 μl of collagenase (1 mg ml^−1) and placed at 37 °C for 15 min. Tissue was pipetted up and down every 2 min, first with a P1000, and then with a P200, until few to no clumps of tissue were visible. For P8 and adult samples, the cortex from an equal number of male and female C57BL/6 mice was dissected and pooled (10 brains for the P8 sample, 4 brains for the adult sample), and cells were dissociated using the neural tissue dissociation kit P (Miltenyi, 130-092-628). For all timepoints, the cell suspension was pelleted (5 min, 800 × g), washed two times with phosphate-buffered saline (PBS) and then resuspended in 180 μl of MACS PEB buffer (PBS (pH 7.2), 0.5% bovine serum albumin and 2 mM EDTA). The cell suspension was then incubated for 15 min with 20 μl of anti-CD31 MicroBeads (Miltenyi Biotec, catalog number 130-097-418) at 4 °C. Cells were then washed with 1 ml of PEB buffer, centrifuged at 300 × g for 5 min and applied to an MS column (Miltenyi, 130-042-201) on a magnetic stand. After three consecutive washes on a magnetic stand with PEB, cells were collected in 0.5 ml of PEB and then pelleted at 300 × g for 10 min at 4 °C. Cells were then resuspended in 1× PBS at a concentration of 50,000 cells per 50 μl, with a viability ≥90% as determined by trypan blue staining, and then used for downstream applications (see below). scRNA-seq of mouse brain cells scRNA-seq libraries were generated using the 10x Chromium Single Cell 3′ v2 reagent kit, according to the manufacturer’s instructions, and were sequenced on an Illumina NextSeq 500. Briefly, raw sequencing data were handled using the 10x Genomics Cell Ranger software ([356]www.10xgenomics.com). Fastq files were mapped to the mm10 genome, and gene counts were quantified using the Cell Ranger count function. Subsequently, expression matrices from each experiment were merged and then imported into Seurat (version 4.0.4)^[357]117 ([358]https://satijalab.org/seurat/) for log normalization. Cells with a percentage of mitochondrial reads above 10%, and with fewer than 250 features, were filtered out. Batch effects were corrected by regressing out the number of mitochondrial read percentage using the vars.to.regress function. Doublet contamination was removed using DoubletFinder (version 2.0.3)^[359]182. Sample integration was achieved using SCTransform (version 0.3.2)^[360]183 before running principal component analysis was performed and significant principal components were used for graph-based clustering. Uniform Manifold Approximation and Projection (UMAP) was used for two-dimensional visualization ([361]https://github.com/lmcinnes/umap). To correct for batch effects, we used Harmony 1.0 ([362]ref. 184). For finding clusters, we used the FindClusters function with a resolution of 0.3. Differential expression of genes per cluster was performed using the Wilcoxon rank sum test (FindMarkers function default). To identify transcription factors regulating the changes in gene expression across endothelial development, we use SCENIC (version 1.2.4)^[363]127. Both ±500 bp and ±10 kb around the mouse TSS motif ranking databases were used for the analysis with default parameters. Cell cluster identities were determined using previously established marker genes validated in the literature^[364]185 or the online tool CellKb ([365]https://www.cellkb.com)^[366]186 leading to the identification of astrocytes or radial glia (Mdk1^+, Igfbp5^+, Meis2^+, Eno1^+)^[367]187, ECs (Cldn5^+, Pecam1^+, Slco1c1^+), pericytes (Abcc9^+, Kcnj8^+, Pdgfrb^+, High1b^+, Vtn^+, Cspg4^+) and SMCs (Acta2^+, Tagln^+)^[368]188, macrophages or megakaryocytes (Pf41^+, Gp1bb^+, Itga2b^+, Clec1b^+), and microglia (Cx3cr1^+, Ccl3^+, Csf1r^+)^[369]189. To identify receptor–ligand interactions, we subset the endothelial, pericyte and microglia clusters from E9.5 and adult mice. The Wilcoxon signed ranked test was used to identify differentially expressed genes between timepoints in each cluster. Only genes present in at least 10% of cells and with a log fold change above 0.25 were considered. We then used CCInx (version 0.5.1)^[370]107 to identify interactions between cell populations across time. Results can be accessed at the interactive shinyapp ([371]https://mcantug.shinyapps.io/Endo_CCInxE9Ad/). Upstream regulation of differentially expressed genes in E9.5 and adult samples was analyzed and visualized by circus plot using NicheNET (version 1.0.0)^[372]92 with default parameter settings. Only active ligands at the 95th quantile were shown. Pseudotime analysis of scRNA-seq For pseudotime analysis, we used Monocle 3 (version 1.0.0)^[373]146. To split the trajectories, we subset the main endothelial Seurat object by cluster. We then converted those objects into Moncole3 compatible cell_data_set (cds) object classes. We processed them according to default parameters. For root node selection, we purposefully chose nodes within the mitotic cluster for each trajectory. For projecting the individual trajectory pseudotime scores calculated in Monocle 3 back onto our original UMAP embedding ([374]Fig. 5), we extracted the pseudotime scored from the Monocle 3 cds and added it to the Seurat object as cell metadata. To calculate the genes whose expression changes across pseudotime, we used the Monocle 3 graph test with a ‘prinicipal_grap’ setting. We subset the gene output that had a q value < 1 × 10^−100. Gene modules were detected with the ‘find_gene_module’ function using default parameters. We then extracted all the genes from each individual gene module and used the AddModuleScore function in Seurat to perform co-expression module scoring on all single cells in the main dataset. Briefly, a gene-expression module score is determined by calculating the average expression of each gene within a list subtracted by the average expression of randomly selected control genes (there is no test statistic, but an aggregated expression score for the input gene list). For gene enrichment analysis, we used Metascape^[375]190 ([376]www.metascape.org). All gene enrichment procedures within metascape were performed with default statistical parameters (minimum overlap > 3 genes, P < 0.01, minimum enrichment > 1.5). Omni-ATAC and RNA-seq of BBB hCMEC/D3 cells Immortalized hCMEC/D3 (Millipore, SCC066) cells were grown to confluence using endothelial cell medium (ScienCell, number 1001) on plates coated with Collagen Type I Rat Tail (Sigma-Aldrich, number C7661). Passages four to six were used for experiments. ATAC libraries were processed as previously described^[377]191. The quality of purified DNA libraries was assessed using the Agilent High Sensitivity DNA kit (Agilent Technologies). Paired-end, 2 × 75 bp sequencing was performed on an Illumina NextSeq 500 instrument. Reads were mapped to the GRCh38 version of the human genome using Bowtie2 with default paired-end settings^[378]170. Mitochondrial reads, reads with a MAPQ < 10 and reads that did not align with the reference genome were removed using Samtools (version 1.13)^[379]171. Duplicated reads were then removed with Picard MarkDuplicates^[380]172. Peak calling was carried out with MACS2 (callpeak –nomodel –broad) (v2.2.7.1)^[381]173. Diffbind (version 3.2)^[382]174,[383]175 was used to import peaksets (min.overlap= 0.66) into RStudio Server (version 1.4.1717, [384]https://www.rstudio.com) using R (version 4.1; [385]ref. 176). The dba.blacklist function was used to filter out peaks that overlap with the ENCODE blacklist. Consensus peaks were converted to mm10 using the LiftOver tool available from the UCSC Genome Browser ([386]https://genome.ucsc.edu/cgi-bin/hgLiftOver). A region was considered conserved if a minimum 0.95 ratio of bases remapped to the mouse genome. Selected regions were also examined using the ECR Browser^[387]192 where the regions were analyzed using rVista 2.0 ([388]ref. 193) to identify conserved transcription factor motifs. The TRANSFAC professional V10.2 vertebrate library was used with default parameters. RNA was isolated using Trizol. Upon processing, RNA from all samples was thawed, and following confirmation of integrity and concentration using a Bioanalyzer, 100 ng was used for low-input library preparation using the NEBNext Ultra II RNA Library Prep kit for Illumina. The libraries were then quantified and sequenced using an Illumina NovaSeq 6000 at a depth of 20 million reads per sample. Reads were first mapped to the human genome (GRCh38) using Salmon (version 1.5.1)^[389]179. Transcript level quantification was then imported using txtimport (version 1.20.0)^[390]180 and analyzed using DESeq2 ([391]ref. 177). Genes were kept and considered actively expressed if they had more than 10 raw counts and >2 log[2](FC) normalized expression. In situ hybridization and confocal imaging Timed matings between C57BL/6 mice were set up, with noon of the day a vaginal plug was detected considered day 0.5. E12.5 embryos (sex undetermined) were collected, and embryonic tissues, as well as postnatal day 7 (P7) mouse brains from female and male mice, were immersion fixed in freshly made 4% PFA, 1× PBS and 0.1% Tween-20 overnight at 4 °C with gentle agitation. Three-week-old and adult brains from male C57BL/6 mice were processed identically, with the exception that they received a retrorbital (100 μl) and transcardiac (150 μl) injection of fluorescently labeled lectin (Vector Laboratories, number DL-1178-1) before tissue collection. After overnight fixation, tissues were washed extensively before being mounted in 1% agarose and sectioned at 35 μm thickness using a Precisionary Instruments Compresstome (VF-300-Z). Sections were then mounted on Superfrost Plus Microscope Slides (Fisher Scientific, number 22-037-246) and then dried at 50 °C for 2 h. Slides were then washed with PBS for 5 min before being allowed to dry for another hour at 60 °C. A hydrophobic barrier around each of the samples was created using the ImmEdge Hydrophobic Barrier Pen (Advanced Cell Diagnostics, catalog number 310018). Sections were then fixed with 10% formalin for 15 min at 4 °C and then underwent successive dehydration steps using 50%, 70% and 100% ethanol. Each dehydration step was performed at room temperature and lasted for 5 min. Following ethanol dehydration, antigen retrieval was then performed using RNA Target Retrieval Reagent (Advanced Cell Diagnostics, catalog number PN 322000). Following antigen retrieval, samples were again serially dehydrated into absolute ethanol. The samples were then incubated in Protease III at 40 °C for 30 min before washing 5 times in ddH[2]O for 1 min each wash. Sections were then incubated in the desired probe (~5–7 drops per section) for 2 h at 40 °C. The slides were then washed 2 times in prepared 1× wash buffer (Advanced Cell Diagnostics, catalog number 310091) for 5 min each wash. Amplification and detection were performed using the RNAscope Multiplex Fluorescent Reagent Kit v2 (Advanced Cell Diagnostics, catalog number 323100). Sections were incubated with Amp1 for 30 min at 40 °C and then washed 2 times in 1× wash buffer (Advanced Cell Diagnostics) for 5 min. This process was repeated for Amp2 and Amp3. Horseradish Peroxidase for Channel 1 (Advanced Cell Diagnostics) was incubated with the sample for 15 min at 40 °C, and then the slides were washed two times in 1× wash buffer for 5 min each wash. Bound probes were detected by adding Opal Dye (Akoya, FP1496001KT) at a 1:1,000 dilution for 30 min at 40 °C, followed by washing two times in wash buffer for 10 min each wash at room temperature. Samples were then incubated in horse radish peroxidase blocking buffer for 30 min at 40 °C followed by two washes in 1× wash buffer for 10 min each wash. Following the RNA Scope steps, samples were then incubated in freshly prepared blocking buffer (10% goat serum, 1× PBS and 0.2% Triton X-100) for an hour at room temperature before adding Alexa Fluor 488-conjugated anti-Claudin-5 (Thermo Fisher; number 4C3C2, number 352588) (1:100) diluted in blocking reagent overnight at 4 °C. The following day, samples were washed extensively in 1× wash buffer (Advanced Cell Diagnostics, catalog number 310091) before being stained with DAPI (ACD Biotechne, catalog number 320858) and mounting with Prolong Gold Antifade Mounting Medium (Thermo Fisher, catalog number [392]P36930). Each assay was processed alongside a positive control slide probed for peptidyl-prolyl isomerase B (ACD Biotechne, catalog number 540651) and a negative control slide treated with dapB (ACD Biotechne, catalog number 707351). Images were collected at a step size of 1 μm for a total stack size of 22 virtual slices using an LSM780 confocal laser scanning microscope (Zeiss) with a ×40 (NA = 1.4) oil objective; laser power = 0.8% (405 nm), 15.7% (488 nm), 2.0% (560 nm) and 5.5% (650 nm); and pinhole = 51 μm. A total of three images were taken per developmental timepoint. Stacks of 20–25 μm, at 0.2 μm per section, were acquired and then compressed to a maximum intensity projection using ZenBlack imaging software (Zeiss). Representative images from each developmental stage were exported as TIFs with individual and merged channels and exported to Adobe Illustrator to generate figures. The number of animals imaged is indicated in the figure legend, with male and female mice being used for the P7 analysis and only male mouse brains for the adult analysis. Lentiviral cloning A Gateway compatible middle entry (pME) plasmid containing luciferase followed by a carboxy terminal P2A cleavage peptide sequence and then a nuclear localized H2A mCherry reporter (JDW 926, pME-Luc-P2A-H2A-mCherry)^[393]194 was used as a template for subcloning each individual transcription factor of interest. Briefly, JDW 926 (pME-Luc-P2A-H2A-mCherry) was digested with BamHI and KpnI to remove luciferase, and then the transcription factor of interest, including a V5 epitope tag followed by a flexible GlyGlyAlaGlyGly (G[2]AG[2]) linker upstream of the open reading frame of the transcription factor, was amplified by PCR, digested to create compatible DNA overhangs with BamHI and KpnI, and then ligated to JDW 926 to create pME-V5-TF-P2A-H2A-mCherry. In the instances where a template was not available, a gene block was directly synthesized (Twist Biosciences) with 20–25-bp homology arms for cold fusion cloning (SBI Biosciences, MC100B-1) into JDW 926. The primers used, templates for PCR or synthesized gene blocks (Twist Biosciences) are listed in [394]Supplementary Table 47. JDW 947/pME-dnTCF4-P2A-mCherry contains a myc-tagged mouse dominant negative TCF4 (TCFL72) lacking the amino terminal 56 amino acids of TCF4, which mediate binding to β-catenin (the nuclear localization sequence and DNA binding domain in the c-terminus of TCF4 are intact, so it can still bind DNA), also referred to as myc-mTcf4ΔN. Of note, there is a missense mutation at the 5′ end, of AAC>GAC, which changes an N to a D residue; otherwise, the entirety of the insert is confirmed. A less conventional MYC tag (MQKLISEEDL) followed by a long linker is present upstream of mouse TCF4. JDW 1391/pME-FLAG-FOXF1-P2A-H2A-mCherry was amplified from Addgene number 153122 (pCS2-FLAG-hFOXF1) (a gift of S. Koch). All TCFs, GATAs, MAFs and FOXF2 middle-entry clones were LR recombined into pLenti6.2-CMV-DEST-PGK-BLAST (JDW 385) (Invitrogen, V36820) to generate lentiviral plasmids for packaging. FOXF1 and ZIC3 were LR recombined into pLenti-CMV-DEST-WPRE-SV40pro-BLAST (Addgene, number 17451) (JDW 1347). All plasmids will be available via Addgene: [395]https://www.addgene.org/Joshua_Wythe/. Lentivirus production and HUVEC transduction For viral packaging, HEK-293T cells (ATCC-CRL3216) were maintained in DMEM, high glucose supplemented with 2 mM GlutaMAX (Gibco, 10566016), 10% FBS (HyClone Characterized FBS, New Zealand origin, 1,000 ml, irradiated 25-40kGy, number SH30406.02) and 1% Pen–Strep (Gibco, 15140122) on 10-cm plates (Corning, 430167). An endotoxin-free 10 μg mixture of equal amounts of the transfer plasmid containing a CMV promoter driving expression of the transcription factor of interest followed by a P2A cleavage sequence and then an H2A-mCherry reporter, along with three third-generation lentiviral packaging plasmids (pMDLg/pRRE containing Gag and Pol, Addgene, number 12251; pRSV-REV, Addgene, number 12253; and pMD2.G encoding VSV-G Env, Addgene, number 12259) was initially diluted in 250 μl of Opti-Mem (Gibco 31985962) before being further diluted in a 1:3 μg ratio with a 1 mg ml^−1 solution of 25 kDa linear polyethylenimine (PEI) (Thermo Fisher, 043896.01). This mixture was then brought to a total volume of 500 μl with Opti-Mem. The DNA-PEI solution was incubated at room temperature for 20 min. Following incubation, the mixture was carefully added dropwise to cells. The plates were gently rocked back and forth to ensure an even distribution of the transfection mixture and then returned to a water-jacketed cell culture incubator set at 37 °C and 5% CO[2]. After 6 h, the medium was aspirated and replaced with antibiotic-free medium (DMEM (Gibco, 10566016) and 10% FBS (HyClone, number SH30406.02). Virus-containing conditioned medium was collected every 24 h post-transfection, over a total period of 72 h, and then centrifuged at 500 × g for 5 min at 4 °C before the supernatant was loaded onto a 10-kDa molecular weight cut-off centrifugal filter unit (Amicon UFC9010) and centrifuged at 500 × g for 35 min at 4 °C to render a final concentration of 20× (from 30 ml to 1.5 ml). The resulting conditioned medium was then aliquoted and snap frozen in liquid nitrogen and stored at −80 °C. Thawed concentrated conditioned medium was then added to low-passage HUVECs (ATCC, CRL-1730) as described below. TEER measurements Each well of a 24-well Transwell plate (Corning 3378) was incubated with 50 μl of rat tail collagen I (Corning CLS354236), diluted to a final concentration of 2 mg ml^−1 in 1× PBS, for 1 h at room temperature. Following incubation with collagen, the collagen was aspirated and the Transwell plates rinsed twice with 1× PBS (Gibco, 10010023). After the Transwell plates were coated with collagen, 500 μl of complete endothelial cell growth medium was added to the lower chamber (Lonza, CC-3162) and 100 μl of the same medium was added to the top chamber to maintain the hydration of the Transwell membrane. Low-passage HUVECs were dissociated using TrypLE (Gibco, 12605010) and then neutralized by the addition of complete endothelial cell medium, after which 100 μl of a cell suspension (at a concentration of 5,000 cells μl^−1, totaling 500,000 cells) was added to the top chamber. HUVECs were maintained on Transwell plates until confluence (defined as TEER of 100 ± 5 Ω cm^−2) was maintained for over 24 h. TEER values were recorded at the same time each day following a 15-min equilibration period at room temperature. After confluence was confirmed, conditioned viral medium was premixed with complete endothelial cell media (ECM) before being added to the top chamber of the Transwell plate. Resistance measurements were taken every 6 h following viral transduction for a total of 48 h. TEER measurements were initially adjusted by subtracting the average values obtained from three Transwell plates that were not seeded with HUVECs. These adjusted measurements were then compared with the non-transduced condition. Each measurement was performed in triplicate per well, with three technical replicates analyzed for each of the three biological replicates. Impedance measurements were taken using the EVOM^2 and the STX3 attachment (World Precision Instruments 91799). FITC–dextran permeability assay To assess cellular permeability, HUVECs were seeded in Transwell plates (Corning 3378) pre-coated with 50 μl of 2 mg ml^−1 collagen I (Corning CLS354236) diluted in 1× PBS, (Gibco, 10010023) as described above. In total, 500,000 low-passage HUVECs were added to the top chamber of the insert and complete endothelial cell growth medium (Lonza, CC-3162) was added to both the upper and lower chambers. Cells were grown for 2 days to create a monolayer and then transduced with lentivirus containing the transcription factors as described above. After 48 h, 10 μl of 1 mg ml^−1 40 kDa fluorescein isothiocyanate (FITC)–dextran (Sigma-Aldrich FD40S) was added to the upper chamber containing the cells and then the cells were incubated for 30 min at 37 °C and 5% CO[2]. Then, 10 μl of medium was collected from the bottom chamber and transferred into a 96-well plate (Corning, CLS3788). Fluorescence (excitation 485 nm, emission 535 nm) was measured using a Cytation TM 1 (BioTek, Agilent). Values are shown as arbitrary fluorescence units. Measurements were initially adjusted by subtracting the fluorescence values of each condition from the fluorescence value produced from the HUVEC medium. These adjusted measurements were then compared with the non-transduced condition. Each measurement was performed in triplicate per well, with three technical replicates analyzed for each of the three biological replicates. Immunofluorescence and confocal imaging of HUVECs HUVECs grown on Transwell membranes for permeability analysis were fixed onto the Transwell plates at the end of the experiment by incubating the inserts in a 24-well dish filled with 4% PFA (pH 7.4) and 1× PBS for 15 min at room temperature. Cells on the Transwell plates were then permeabilized in 0.1% Triton X-100 and 1× PBS for 30 min at room temperature, rinsed twice in 1× PBS at room temperature and then incubated in blocking buffer (1% bovine serum albumin (Fraction V), 0.1% Triton X-100 and 1× PBS) for 1 h at room temperature. The wells were then incubated with primary antibodies diluted in blocking buffer (at the concentrations indicated below) in a humidified chamber overnight at 4 °C. The following day, the cells were washed three times, 5 min per wash, in 1× PBS with gentle agitation on an orbital shaker at room temperature, and then incubated in Alexa dye-conjugated secondary antibody diluted in blocking buffer for 1 h at room temperature. The inserts were then washed for 5 min, three times, in 1× PBS, and then incubated in DAPI (diluted 1:1,000 in 1× PBS) to label nuclei. Finally, the inserts were washed for 5 min, three times, in 1× PBS. Transwell inserts were inverted and removed with a scalpel and then mounted on a glass slide (Fisher, 22-037-246) in Prolong Gold Antifade mounting medium (Thermo Fisher, [396]P36930) and coverslipped. Antibodies and dyes used are listed in [397]Supplementary Table 47. Afterward, a total of 3 random 2,480 × 2,480 μm per pixel images per condition were obtained from 3 biological replicates on a ZEISS LSM 700 laser scanning microscope with a ×40 (NA = 1.4) oil objective; laser power = 2.0% (405 nm), 18.0% (488 nm), 8.0% (560 nm) and 2.5% (650 nm); and pinhole = 51 μm. Z-stacks of 10–15 μm, 0.2 μm per section, were collected and then compressed into a maximum intensity projection using ZenBlue imaging software (Zeiss). Representative images from each condition were then exported as both individual and merged channels as.TIF files to Adobe Illustrator. Representative images are shown. qRT–PCR A total of 3.0 × 10^6 low-passage HUVECs were seeded onto collagen I-coated 10-cm plates (Corning, 430293) and cultured for 2 days in complete ECM medium. At 60% confluence, concentrated lentivirus-conditioned medium from HEK-293T cells was added to the HUVECs at a volumetric ratio of 1:5 conditioned medium to complete ECM. After 2 days, when nuclear mCherry expression was evident, cells were rinsed twice in 1× DPBS (Gibco, 14190144), and then lysed in 1,000 μl of Trizol (Invitrogen, 10-296-010). Lysates were collected and RNA extracted using the RNeasy Plus Kit (Qiagen, 74104). Total RNA was quantified by Nanodrop (Thermo Scientific ND-2000C) and 1 μg of total RNA was reverse transcribed in a 20-μl reaction using the high-capacity cDNA reverse transcription kit (Applied Biosystems, 43-688-14). The cDNA was then diluted 12.5-fold and 50 ng of the cDNA was subsequently used in a 20-μl qPCR with Sybr Green PCR Master Mix (Applied Biosystems 4309155). Data were normalized to TATA-box binding protein (TBP) using the ΔΔCt method^[398]195 and are shown as fold change (calculated as 2^−ΔΔCt). Technical triplicates were analyzed for each biological sample, with three biological replicates (from different transduction experiments) analyzed per condition. Data were analyzed using PRISM software (GraphPad) and shown as the mean ± s.d. Statistical significance was determined using two-tailed Student’s t-test (for two groups) or a multi-way ANOVA (for multiple comparisons with more than two groups) with P < 0.05 considered significant. Primers used for real-time quantitative reverse transcriptase polymerase chain reaction (qRT–PCR) were designed with the following criteria: ~100–250 bp size amplicon, span at least one exon–exon junction (if possible), generate a single amplicon (as determined by melt curve analysis and agarose electrophoresis) and show an efficiency above 90%. Primers are listed in [399]Supplementary Table 49. General statistics Unless otherwise indicated, experiments were performed using a minimum of two independent biological replicates. The exact number of animals used per biological sample is indicated in each experimental method described above (for example, ‘Mouse endothelial nuclear isolations’ and the processing sections for the RNA-seq or ATAC-seq descriptions). Calculations to determine statistical significance are detailed within each experimental method section. For embryonic and adult mouse samples, a mixture of female and male samples was collected for all analyses. Reporting summary Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article. Extended Data Extended Data Fig. 1 |. Lymphatic marker expression in RNA-Seq and ATAC-Seq. Extended Data Fig. 1 | [400]Open in a new tab a) Expression levels of endothelial markers (top) and lymphatic markers (bottom) as determined by RNA-Sequencing of adult endothelial and input samples (normalized expression). b) Genome browser tracks of lymphatic markers in all ATAC-seq tissue samples reveals no substantial difference in genome accessibility at loci for lymphatic markers between input and endothelial enriched samples. Extended Data Fig. 2 |. Transcription factor expression and motif presence across all organs. Extended Data Fig. 2 | [401]Open in a new tab a) Expression levels of common transcription factors (related to [402]Fig. 2) across the adult endothelium (normalized expression). b) Enriched motifs identified by HOMER (one-tailed Fisher’s exact test) from all organs, with all timepoints condensed into one sample per organ. Size of the bubble and the color represent the p-value. The top 50 motifs are shown (values represent the −log[10] of the p-value). c) Expression levels of enriched transcription factors across the adult endothelium (normalized expression). Extended Data Fig. 3 |. Chromatin Accessibility Changes Across Time in the Heart Endothelium. Extended Data Fig. 3 | [403]Open in a new tab a) Differential chromatin accessibility determined by ATAC-Seq peaks in the heart endothelium (11,079 peaks) at E12.5, postnatal day 6 (P6) and adult (2-month-old) mice. b) Biological processes from expressed genes and with accessible chromatin in each timepoint. c) Top 10 transcription factor motifs ranked by gene expression for each age. Log2 expression over input indicated in the y-axis. Motif enrichment p-value (one-tailed Fisher’s exact test) is shown according to the dot size. Extended Data Fig. 4 |. Chromatin Accessibility Changes Across Time in the Liver Endothelium. Extended Data Fig. 4 | [404]Open in a new tab a) Differential chromatin accessibility determined by ATAC-Seq peaks in the liver endothelium (8,666 peaks) at E12.5, postnatal day 6 (P6) and adult (2-month-old) mice. b) Biological processes from expressed genes and with accessible chromatin in each timepoint. c) Top 10 transcription factor motifs ranked by gene expression for each age. Log2 expression over input indicated in the y-axis. Motif enrichment p-value (one-tailed Fisher’s exact test) is shown according to the dot size. Extended Data Fig. 5 |. Chromatin Accessibility Changes Across Time in the Lung Endothelium. Extended Data Fig. 5 | [405]Open in a new tab a) Differential chromatin accessibility determined by ATAC-Seq peaks in the lung endothelium (1,731 peaks) at E12.5, postnatal day 6 (P6) and adult (2-month-old) mice. b) Biological processes from expressed genes and with accessible chromatin in each timepoint. c) Top 10 transcription factor motifs ranked by gene expression for each age. Log2 expression over input indicated in the y-axis. Motif enrichment p-value (one-tailed Fisher’s exact test) is shown according to the dot size. Extended Data Fig. 6 |. Chromatin Accessibility Changes Across Time in the Kidney Endothelium. Extended Data Fig. 6 | [406]Open in a new tab a) Differential chromatin accessibility determined by ATAC-Seq peaks in the kidney endothelium (3,035 peaks) at E12.5, postnatal day 6 (P6) and adult (2-month-old) mice. b) Biological processes from expressed genes and with accessible chromatin in each timepoint. c) Top 10 transcription factor motifs ranked by gene expression for each age. Log2 expression over input indicated in the y-axis. Motif enrichment p-value (one-tailed Fisher’s exact test) is shown according to the dot size. Extended Data Fig. 7 |. Novel potential gene regulatory elements in organ specific transcription factors. Extended Data Fig. 7 | [407]Open in a new tab Genome browser tracks from ATAC-seq data highlight transcription factor DNA-binding motifs in accessible chromatin regions enriched in endothelial samples from the mouse (a) brain, (b) heart and (c) liver. Motifs present and their location are highlighted in black, while the genomic loci are indicated below each panel in blue (boxes represent exons, arrows indicate the direction of transcription). Extended Data Fig. 8 |. Cell to Cell Communication Changes in the Neurovascular Unit Over Time. Extended Data Fig. 8 | [408]Open in a new tab a) Schematic representation of the harvesting and isolation of endothelial cells from E9.5, E12.5, E16.5, P8 and adult mice. Cells were purified using Magnetic Isolation Cells Sorting (MACS) and processed for downstream sequencing and analysis following the 10x Genomics protocol. b) UMAP of cells clustered by collection timepoint. c) UMAP of cells with the identity of each cluster indicated. d) Feature plot showing expression of the endothelial marker Cdh5 (encodes VE-Cadherin) superimposed on the UMAP to identity endothelial cells (normalized expression values). e) Circos plot of differentially expressed ligands in non-EC cells within our dataset, as well as their target genes expressed in the CNS endothelium between E9.5 and Adult. f) Unbiased analysis of top predicted interactions of differentially expressed ligands and receptors between ECs and pericytes in E9.5 and adult using the Cell-Cell Interactions (CCInx using adjusted p-value, and average log2 fold change). Wilcoxon rank sum test, two-sided P-value. Extended Data Fig. 9 |. Differential gene expression in endothelial cells across time. Extended Data Fig. 9 | [409]Open in a new tab a) UMAP of only endothelial cells (from [410]Fig. 6) colored coded by timepoint and clustered using Monocle 3 to show maturation across time. b) Feature plots denoting normalized expression of endothelial marker genes, Pecam1 and Cdh5, and canonical lymphatic marker genes, Lyve1 and Prox1, demonstrate the preferential enrichment of blood endothelial cell transcriptomes rather than lymphatic endothelial transcriptomes. c) Heatmap of differential gene expression analysis of single cell RNA-seq results from brain endothelial cells profiled across developmental time (row z-score). Genes in red have a known role in blood brain barrier function. The top 10 differentially expressed genes for each time point are shown, followed by the pan-endothelial transcripts Pecam1, Cldn5 and Kdr (normalized expression). d) Expression of BBB-related genes over developmental time. Extended Data Fig. 10 |. Gene Ontology analysis of endothelial maturation within the mouse brain. Extended Data Fig. 10 | [411]Open in a new tab The differentially expressed gene signatures from the single cell endothelial subtypes of the developing and adult mouse brain that were used for pseudotime analysis by Moncole3 were processed for Gene Ontology analysis. The heatmap reflects which pathways and processes are upregulated in each of the unique endothelial cell subtypes (values represent the −log[10] of the p-value). Colored by significance of pathway enrichment (2-sided p-value). Supplementary Material Cantu_Gutierrez_et_al_Supplemental_Information_v3 [412]NIHMS2069201-supplement-Cantu_Gutierrez_et_al_Supplemental_Informa tion_v3.pdf^ (1.5MB, pdf) Cantu_et_al_Supplementary_Tables_v2 [413]NIHMS2069201-supplement-Cantu_et_al_Supplementary_Tables_v2.xlsx^ (12.6MB, xlsx) Acknowledgements