Blood vessel growth and remodelling are essential during embryonic development and disease pathogenesis. The diversity of endothelial cells (ECs) is transcriptionally evident and ECs undergo dynamic changes in gene expression during vessel growth and remodelling. Here, we investigated the role of the histone acetyltransferase HBO1 (KAT7), which is important for activating genes during development and for histone H3 lysine 14 acetylation (H3K14ac). Loss of HBO1 and H3K14ac impaired developmental sprouting angiogenesis and reduced pathological EC overgrowth in the retinal endothelium. Single-cell RNA sequencing of retinal ECs revealed an increased abundance of tip cells in Hbo1-deficient retinas, which led to EC overcrowding in the retinal sprouting front and prevented efficient tip cell migration. We found that H3K14ac was highly abundant in the endothelial genome in both intra- and intergenic regions, suggesting that HBO1 acts as a genome organiser that promotes efficient tip cell behaviour necessary for sprouting angiogenesis.
The expansion of blood vessel networks by angiogenesis is crucial for matching blood supply to the metabolic demands of growing tissues. In adult tissues, blood vessels and the endothelial cells (ECs) that line them generally remain in a quiescent state. However, the reactivation of angiogenesis is a hallmark of a number of diseases, including cancer and vision-threatening eye diseases (Ferrara and Adamis, 2016; Miller et al., 2013). Sprouting angiogenesis and subsequent vessel remodelling generates diverse EC identities (Potente and Makinen, 2017) and is a major driver of vessel growth in both developmental and pathological contexts. Pro-angiogenic signals such as vascular endothelial growth factor (VEGFA) stimulate the formation of vessel sprouts. Vessel sprouts are made up of two types of EC: migratory tip cells that guide growing sprouts and trailing stalk cells that proliferate to contribute new cells to growing sprouts (Gerhardt et al., 2003). Tip and stalk cells exist transiently, switching identities during angiogenesis (Jakobsson et al., 2010; Arima et al., 2011), with tip cells ultimately directed to the growing arteries (Pitulescu et al., 2017; Xu et al., 2014). The new vessels formed during sprouting angiogenesis are immature and must be remodelled into hierarchical, mature vessel networks. This includes establishment of arteries, veins and capillaries and defined identity of ECs within each of these (Potente et al., 2011).
The diversity of these distinct EC types is evident at a transcriptional level, with each type exhibiting gene expression profiles that distinguish them from another. Large-scale characterisation of transcriptional heterogeneity of EC populations has been made possible by advances in single-cell RNA sequencing (scRNA-seq) (Jakab and Augustin, 2020). This has enabled profiling of EC heterogeneity in a range of tissues, including within quiescent adult tissues (Vanlandewijck et al., 2018; Kalucka et al., 2020; Tikhonova et al., 2019; Lukowski et al., 2019), in the activated endothelium of tumours (Goveia et al., 2020; Zhao et al., 2018) and in some developing vascular beds, including in the brain (Sabbagh et al., 2018) and coronary arteries (Su et al., 2018). These studies have also highlighted key transcription factors that are expressed within each subtype. The actions of DNA-binding transcription factors are essential for cell type-specific gene expression; however, the compacted nature of chromatin is refractory to these interactions. This requires the additional actions of chromatin-associating enzymes that modify chromatin. Snapshots of EC heterogeneity obtained from scRNA-seq do not provide information about the dynamic interactions of these enzymes with chromatin that enable gene expression patterns to be established. Furthermore, the plasticity and fate switching that occurs in the immature vessel networks to enable efficient vessel network expansion and remodelling requires that ECs change their gene expression patterns. This is particularly relevant in the context of tip and stalk cell shuffling, as entering or exiting the transcriptionally distinct tip cell fate would involve dynamic changes in gene expression (del Toro et al., 2010; Strasser et al., 2010).
Histone acetyltransferases (HATs) are chromatin-modifying enzymes that dynamically regulate gene expression by catalysing the transfer of an acetyl group to histone tail lysine residues. The resulting histone acetylation mark is typically correlated with active gene expression (Karmodiya et al., 2012; Wang et al., 2009). Various histone lysine residues can be acetylated with differences in genome abundance, function and downstream consequences. Here, we focussed on HBO1 (also known as KAT7 or MYST2), which is exclusively required for histone H3 lysine 14 acetylation (H3K14ac) (Kueh et al., 2020, 2011), the second most abundant mark in the human genome with >15% of sites acetylated (Feller et al., 2015; Hansen et al., 2019). Of interest, HBO1 has an essential role in regulating the gene expression patterning necessary for determining embryonic structures and cell types during development (Kueh et al., 2011). Deletion of Hbo1 leads to a global loss of H3K14ac and a failure of de novo expression of lineage-specific genes required for post-gastrulation differentiation. This includes failure to express the genes Tie1, Tie2 (Tek) and Vegfr2 (Kdr), which are necessary for patterning the vasculature (Kueh et al., 2011). In addition, Hbo1–/– embryos display a range of defects, including but not limited to abnormal yolk sac vasculature and dilated dorsal aortae and pharyngeal arch arteries (Kueh et al., 2011). In zebrafish, the Hbo1 orthologue kat7b was shown to be required for normal expression of Vegfr2 and necessary to promote EC migration and vessel sprouting (Yan et al., 2018). Together, the observation that HBO1 is necessary for lineage-specific gene expression accompanying differentiation and the vascular phenotype observed in kat7b morphants and Hbo1–/– embryos posit a role for HBO1 in regulating gene expression patterns required for EC differentiation during angiogenesis.
Using tissue-specific conditional deletion of Hbo1 in mice, scRNA-seq, chromatin immunoprecipitation and sequencing, whole-mount immunofluorescence and an oxygen-induced retinopathy model, we have investigated the role of HBO1 in establishing the identity of ECs involved in angiogenesis in the mouse retina during postnatal development and in a disease model.
HBO1 promotes sprouting angiogenesis
The neonatal mouse retina is a widely used model of angiogenesis, in which both vessel sprouting in response to pro-angiogenic factors and subsequent vessel network remodelling occur in parallel (Fig. 1A). EC subtypes involved in both aspects can be readily identified. To study the EC-autonomous role for HBO1 during sprouting angiogenesis, mice carrying Hbo1 conditional alleles (Kueh et al., 2011) were crossed to those expressing the tamoxifen-inducible, endothelial-specific Cdh5(PAC)-creERT2 transgene (Wang et al., 2010). Pups were administered tamoxifen at postnatal day (P)1 and P2 to induce nuclear translocation of cre-recombinase and EC-specific deletion of exon 1 of Hbo1 (Fig. 1B), hereafter referred to as Hbo1iEC/iEC mice. The genotypes of control mice included Hbo1fl/+;Cdh5(PAC)-creERT2+ (heterozygote following tamoxifen administration), Hbo1fl/+ and Hbo1fl/fl. Deletion of Hbo1 resulted in a 63% reduction in endothelial Hbo1 expression in lung ECs (Fig. 1C). H3K14 is the major acetylation target of HBO1 during mouse embryonic development and in human cell lines (Kueh et al., 2011, 2020). Congruently, Hbo1iEC/iEC mice showed a significant reduction in endothelial H3K14ac compared with control mice (Fig. 1D,E). The P6 retina of Hbo1iEC/iEC mice displayed impaired angiogenic vessel growth as shown by reduced vessel area and radial vessel outgrowth compared with control mice (Fig. 1F-H). Furthermore, the density of the vessel network formed was reduced in Hbo1iEC/iEC mice as shown by fewer vessel branchpoints and vessel segments (Fig. 1I,J). The patterning of arteries and veins was normal in Hbo1iEC/iEC retinas based on known morphological differences (Crist et al., 2017) and we did not observe any incidence of either vessel crossing over the other (Fig. 1F).
Reduced pathological angiogenesis in the absence of HBO1
Similar mechanisms can contribute to both developmental and pathological vessel growth (Dubrac et al., 2016; Raimondi et al., 2014; Rama et al., 2015). We used the oxygen-induced retinopathy (OIR) model to determine whether HBO1 is also important for pathological vessel growth. Mice were exposed to high oxygen for 5 days from P7, leading to vaso-obliteration of the central retina, which induces VEGFA-dependent pathological neovascular growth (Lee et al., 2013; Aiello et al., 1995). The conditional Hbo1iEC/iEC mutation was induced following exposure to high oxygen and retinas were then examined at the peak of pathological neovascular growth at P17 (Fig. 2A). We have previously established that neovascular lesions stain brightly for the basement membrane protein collagen IV (Grant et al., 2020). Staining for collagen IV revealed a significant reduction in pathological neovascular area in Hbo1iEC/iEC mice compared with controls (Fig. 2B,C), suggesting that HBO1 is also required for pathological vessel growth. Consistent with the developmental angiogenesis defect, Hbo1iEC/iEC retinas had a significantly increased vaso-obliterated area compared with controls (Fig. 2B,D), a region that is revascularised by sprouting angiogenesis from the unaffected peripheral vessels. The results from the OIR model support a necessary role for HBO1 in regulating both developmental and pathological vessel growth.
HBO1 is not essential for EC proliferation, EC apoptosis, vessel remodelling or maintaining vessel patterning
HBO1 has been shown to regulate proliferation in other cell types (Johmura et al., 2008; Wu and Liu, 2008). To assess proliferation, P6 mice were labelled for 2 h with the thymidine nucleoside analogue 5-ethynyl-2′-deoxyuridine (EdU), which is incorporated into newly synthesised DNA during S phase of the cell cycle. Proliferation rates in Hbo1iEC/iEC retinas were not significantly different from control mice (Fig. S1A,B). Further supporting this, the number of ECs relative to vessel area across the retina was unchanged between control and Hbo1iEC/iEC mice (Fig. S1C). During sprouting angiogenesis, the initially dense vessel plexus remodels into a hierarchical network involving both apoptosis-dependent and -independent vessel regression (Watson et al., 2017; Grant et al., 2020). We assessed apoptosis by staining for cleaved (active) caspase 3 and found that Hbo1iEC/iEC mice showed normal levels of endothelial apoptosis (Fig. S1D,E). We also found that vessel regression was unchanged in Hbo1iEC/iEC retinas (Fig. S1D,F). These results suggest that defects in proliferation or aberrant vessel pruning do not account for reduced vascularisation in Hbo1iEC/iEC retinas.
HBO1 was not required for vessel maintenance in adult tissues, as deletion of Hbo1 at 6 weeks of age had no effect on the vasculature (Fig. S2A-D). Interestingly, H3K14ac expression was barely detectable in adult ECs from control retinas compared with surrounding cells and this was not reduced further in Hbo1iEC/iEC retinas (Fig. S2E). We confirmed that cre-recombinase was active in these mice by intercrossing the cre-inducible, cell membrane-targeted EGFP reporter allele mTmG (Muzumdar et al., 2007) into these mice, which leads to EGFP expression in ECs, and found extensive EGFP expression in ECs of both control and Hbo1iEC/iEC genotypes (Fig. S2F). Despite reduced H3K14ac in adults, we found that Hbo1 mRNA levels were consistent across all ages from P6-50 in a publicly available developmental retinal EC RNA-seq dataset (Jeong et al., 2017) (Fig. S2G). Taken together, these findings suggest that HBO1 is not required during the neonatal remodelling process or for maintaining vessel patterning in established vessel networks.
Profiling retinal EC transcriptomes by single-cell RNA sequencing
To investigate our hypothesis that HBO1 is necessary for specifying EC identity and to determine whether this could explain the phenotype caused by EC-specific loss of HBO1 in the retina, we first established the transcriptional heterogeneity of ECs within the postnatal retina at the single-cell level. We sorted ECs from dissociated P6 retinas into 384-well plates and performed scRNA-seq using CEL-Seq2 (Hashimshony et al., 2016). For this, we used heterozygous control and Hbo1iEC/iEC mice crossed with the mTmG reporter and treated with tamoxifen at P1 and P2 (hereafter ControliEC/+;mTmG and Hbo1iEC/iEC;mTmG, respectively). We defined ECs as PECAM1+EGFP+, thereby isolating only those ECs in which cre recombination had taken place. We sequenced 1733 cells with an average sequencing depth of 13,222 reads per cell and 3466 genes per cell. After quality controlling for the number of detected genes and mitochondrial and spike-in RNA reads (Fig. S3A), 1393 cells were pooled, batch-corrected, clustered and visualised using uniform manifold approximation and projection (UMAP) plots (Fig. 3A, Fig. S3B,C). All cells expressed high levels of the endothelial markers Pecam1, Cdh5 and Vegfr2 (Fig. S3D). Consistent with the neural nature of the retina, retinal ECs were transcriptionally reminiscent of adult brain ECs, including robust expression of the brain EC-specific marker Pglyrp1 (Kalucka et al., 2020) (Fig. S3E,F). Clustering analysis identified seven clusters (Fig. 3A, Fig. S4A,B) that were annotated based on known EC-subtype markers such as Gja4 for arterial ECs (Fang et al., 2017), Aplnr for venous ECs (Saint-Geniez et al., 2003) and Esm1 for tip cells (Rocha et al., 2014), in addition to recent other endothelial subtype marker gene lists (Kalucka et al., 2020; Goveia et al., 2020) (Fig. 3B, Fig. S4A,B). Capillary ECs exhibit tissue-specific gene expression (Kalucka et al., 2020). Therefore, we annotated capillary ECs using brain capillary markers, such as Spock2 for all capillary ECs (Kalucka et al., 2020), the arterial capillary marker Tgfb2 and the venous capillary marker Tfrc (Vanlandewijck et al., 2018) (Fig. 3B). Three clusters expressed venous/capillary markers and no proliferation markers (Fig. 3A,B, Fig. S4A). For clarity of presentation and downstream analysis, we combined these three clusters. We were unable to discern a discrete stalk cell population even using known stalk cell-expressed genes, including those activated by Notch and TGFβ/BMP signalling (Moya et al., 2012; Larrivee et al., 2012; Strasser et al., 2010), and other genes expressed by stalk-like cells identified by scRNA-seq of tumour ECs (Zhao et al., 2018) (Fig. S4A). In contrast to adult tissue vasculature in which very few ECs were found to be proliferating (Kalucka et al., 2020), we identified two proliferative clusters based on cyclin gene expression that accounted for about 30% of all cells analysed (Fig. S4A-C). Most of the ECs within the proliferative clusters were of venous/capillary identity (Fig. 3B). This is consistent with these ECs being highly proliferative in the postnatal retina (Ehling et al., 2013). Supporting the annotation of distinct EC types, heatmap analysis of top-ranking marker genes revealed distinct transcriptional signatures of ECs within each cluster (Fig. 3C, Table S1). Thus, we were able to interrogate intra-tissue endothelial heterogeneity in the postnatal retina, identifying tip cells, which predominate in the sprouting region, and arterial and venous ECs, which predominate in the remodelling region of the neonatal retina.
Hbo1-deficient retinas have increased expression of tip cell genes and an expanded tip cell compartment
H3K14ac is a mark correlated with active gene transcription (Karmodiya et al., 2012). As HBO1 is essential for H3K14ac in mouse cells (Kueh et al., 2011; Mishima et al., 2011), human cells lines (Kueh et al., 2020) and mouse retinal ECs (Fig. 1D,E), we predicted that loss of HBO1 could affect gene expression. Given the relatively small number of cells available for scRNA-seq per animal in the neonatal retina, which was further reduced when the cells were divided up into clusters, we examined the effects of loss of HBO1 on gene expression ignoring the cluster identity. Neonatal retinal ECs from Hbo1iEC/iEC mice had 92 downregulated and 61 upregulated genes with false discovery rate (FDR) <0.05 compared with ControliEC/+;mTmG ECs (Fig. 4A,B, Table S2). Differentially expressed genes in the venous/capillary cluster (which contained the most cells), accounted for about 40% of the differentially expressed genes identified ignoring cluster identity (Fig. S5A). Appropriate deletion of the Hbo1 gene was confirmed by a reduction in Hbo1 mRNA, and there was no change in the expression of other HAT genes in the absence of HBO1 (Fig. S5B). Mirroring the normal artery and vein patterning in the Hbo1iEC/iEC retina, we found the mRNA levels of artery and vein markers were not consistently up- or downregulated in Hbo1iEC/iEC;mTmG retinal ECs (Fig. S5C,D). In contrast, we found that many tip cell-expressed genes, including Esm1, Igfbp3 and Cxcr4 (del Toro et al., 2010; Strasser et al., 2010), were more highly expressed in Hbo1iEC/iEC;mTmG ECs compared with ControliEC/+;mTmG ECs (Fig. 4C). In addition, Hbo1iEC/iEC;mTmG retinas tended to have a higher proportion of tip cells expressing tip cell markers highly (Fig. 4D, Fig. S5E). We found that upregulated genes were associated with gene ontology (GO) terms for ‘angiogenesis’, ‘blood vessel/vasculature development’ and ‘cell adhesion’ were amongst the most significantly enriched in Hbo1iEC/iEC;mTmG mice (Fig. 4E,F, Table S3). In contrast, many of the downregulated genes encoded transport proteins from solute carrier SLC families (Table S2) and we saw a corresponding predominance of downregulated GO terms relating to transport of small molecules and metabolites (Fig. 4E, Fig. S5F, Table S3). Overall, these data provide insight into the processes that are dependent on HBO1 during angiogenesis. Deletion of Hbo1 did not affect genes in the Notch, VEGF, FOXO, TGFβ and HIF1 signalling pathways in a directional manner (Fig. S5G). To determine whether there was an accumulation or depletion of a particular subtype of EC, we also analysed the differential abundance per genotype within each cluster identified by scRNA-seq. We found that Hbo1iEC/iEC;mTmG retinas had an enrichment of cells in the tip cell cluster (Fig. 4G). Consistent with our EdU analysis showing no proliferation defect in Hbo1iEC/iEC retinas, we found that ECs in each proliferative cluster were similarly abundant in both genotypes (Fig. 4G). Overall, differential analysis of our scRNA-seq data suggested that tip cell genes were upregulated per cell and that there was an accumulation of tip cells in retinas in which Hbo1 had been deleted.
H3K14ac is widely distributed across the endothelial genome but enriched at genes required for vascular development
Previous studies have shown that H3K14ac is a highly abundant mark in the genome (Feller et al., 2015; Hansen et al., 2019). To understand whether its occupancy correlates with particular regions of the activated endothelial genome, we performed chromatin-immunoprecipitation and sequencing (ChIP-seq) in human umbilical vein ECs (HUVECs). To reproduce angiogenic conditions, HUVECs were grown in low oxygen (3% oxygen) and stimulated with VEGFA prior to harvesting for ChIP (Fig. S6A). There were a large number of genes for which H3K14ac was enriched compared with histone H3 (Table S4), including endothelial genes such as PECAM1 and DLL4 (Fig. 5A, Fig. S6B). Supporting this, among the most significant GO biological processes terms for gene bodies enriched with H3K14ac above histone H3 were ‘angiogenesis’, ‘blood vessel development’ and ‘blood vessel morphogenesis’ (Fig. 5B, Table S5), highlighting a crucial role for H3K14ac in behaviours essential for vascular development. H3K14ac-depleted regions were associated with GO terms relating to neural processes, such as ‘regulation of postsynaptic membrane potential’ and ‘GABAergic synaptic transmission’, likely reflecting the repressed state of genes specific to other tissues (Fig. S6C, Table S5). We compared genes with enriched H3K14ac with genes previously found to be highly expressed in activated HUVECs (Zhang et al., 2013). We found that genes highly expressed in HUVECs were significantly enriched for H3K14ac, whereas genes expressed at a low level were depleted for H3K14ac (Fig. 5C). Consistent with the expression pattern of HBO1 occupancy (Saksouk et al., 2009; Avvakumov et al., 2012; Lalonde et al., 2013), we found that the level of H3K14ac immediately downstream of the transcription start site (TSS) correlated positively with gene expression (Fig. 5D). Interestingly, H3K14ac was widely distributed across the genome with only ∼35% reads mapping to gene bodies (Fig. 5E, Fig. S6D).
HBO1 is required for normal tip cell sprouting
As the scRNA-seq suggested an increase in the proportion of tip cells in the Hbo1iEC/iEC retinas (Fig. 4G), we examined tip cell numbers and function in the Hbo1iEC/iEC retinas. Whereas control mice had many vessel sprouts indicative of tip cell activity extending out from the vessel front, Hbo1iEC/iEC retinas had a blunted appearance (Fig. 6A) and a significant decrease in the number of vessel sprouts (Fig. 6B). However, the sprouts that did form were on average the same length in control and Hbo1iEC/iEC mice (Fig. 6C). Following vascularisation of the superficial layer of the retina, vessels sprout downwards to form additional vessel layers deeper within the retina, beginning around P7. Although the superficial layer of Hbo1iEC/iEC retinas was vascularised to the same extent as control retinas by P10 (Fig. S7A,B), sprouting into the deeper layers of the retina was impaired as evidenced by a reduction in the formation of the deep vessel layer in P10 Hbo1iEC/iEC mice (Fig. S7C). These findings confirm a general vessel sprouting defect.
To determine whether tip cell identity was correctly specified in the absence of HBO1, retinas were stained for the VEGFA-responsive marker ESM1, which is selectively expressed by tip cells in the sprouting retinal vasculature (Rocha et al., 2014). We found that ECs at the sprout front expressed ESM1, indicating that they were able to respond to VEGFA signalling and adopt tip cell identity (Fig. 6A). Consistent with this, staining for another tip cell marker, DLL4, showed ECs at the sprouting front of Hbo1iEC/iEC retinas expressing DLL4, consistent with successful tip cell specification (Fig. S7D). ESM1+ ECs were normally distributed to the sprouting front of Hbo1iEC/iEC retinas, but there were significantly more of them than in control mice (Fig. 6A,D). Consistent with the increased Esm1 mRNA levels detected by scRNA-seq, we found that expression of ESM1 protein in Hbo1iEC/iEC ESM1+ cells was increased compared with control ESM1+ cells (Fig. 6E). In contrast, the number of ECs expressing another tip cell marker, hypoxia-responsive ANG2, was not significantly changed across the Hbo1iEC/iEC retina compared with controls (Fig. S7E,F). These results, together with the scRNA-seq data, suggest that the altered sprouting in Hbo1iEC/iEC mice was not due to a failure to specify tip cell identity, but may be due to an inability to form functional sprouts.
Hbo1-deficient tip cells fail to establish biased migration and accumulate in the sprouting front
HBO1 has been previously shown to promote EC migration using in vitro scratch-wound assays (Yan et al., 2018). Tip cells are highly migratory and a defect in this migratory behaviour could explain the reduced vessel sprouting observed in Hbo1iEC/iEC mutants. Tip cells achieve directed migration by establishing front-rear polarity, which occurs as the tip cell orients its Golgi apparatus towards the sprouting front relative to the nucleus (Moya et al., 2012; Franco et al., 2015; Dubrac et al., 2016) (Fig. 6F). To assess Golgi apparatus orientation in tip cells, we stained retinas with the Golgi apparatus marker GM130, the nuclear marker FLI1 and the tip cell marker ESM1. Consistent with previous findings (Dubrac et al., 2016), in control mice, most tip cells were preferentially polarised in the direction of sprout elongation (Fig. 6G,H). In contrast, Hbo1iEC/iEC tip cells did not display a polarity bias and Golgi apparatus were randomly aligned relative to the direction of sprout elongation (Fig. 6G,H), suggesting an impairment in sprout migration. Consistent with failed sprout migration, the sprouting front of Hbo1iEC/iEC retinas displayed increased EC density (Fig. 6G,I), corresponding to an increase in vessel width compared with control vessels (Fig. 6J). This is in contrast to vessels in the remodelling zone, which displayed normal EC density and vessel width (Fig. 6K,L), further supporting normal vessel remodelling in the absence of HBO1. Together, these data suggest that in the absence of HBO1 tip cells fail to sprout and migrate normally, which leads to an accumulation of ECs within the sprouting front and a relative greater abundance of tip cells, as observed in the scRNA-seq data.
Establishing distinct EC identity and function is crucial for the normal growth and patterning of the vasculature and involves changes in gene expression patterns. Although it is known that ECs exhibit a distinct identity at the transcriptional level, how these differences arise through the actions of chromatin-associating enzymes is not well understood. In this study, we investigated the involvement of the histone acetyltransferase HBO1 in EC identity specification. HBO1 is responsible for deposition of the second most abundant histone acetylation mark in the human genome, H3K14ac (Kueh et al., 2011, 2020; Feller et al., 2015; Hansen et al., 2019). We found that H3K14ac was widely distributed across the human endothelial genome, and that deletion of Hbo1 in mice reduced global H3K14ac. In the absence of H3K14ac and HBO1, retinal angiogenesis was impaired in both normal and pathological contexts. This was not due to a failure of tip cell specification but rather to an accumulation of tip cells, which were unable to become organised collectively and efficiently in order to undergo directed migration.
HBO1 was previously implicated in the regulation of angiogenic EC migration in HUVECs and zebrafish (Yan et al., 2018). Our study represents the first to directly examine the role of a HAT during in vivo mammalian vascular development. Other studies have shown that the HAT activity of EP300 is necessary for transcription factor binding and activation of essential endothelial genes (Sacilotto et al., 2016; Shu et al., 2015); however, little is known about the cell-autonomous roles of other HATs during sprouting angiogenesis or vascular development more generally. Loss of the histone deacetylase activity of HDAC10 changes the ability of ECs to sprout in vitro (Duan et al., 2017), supporting the role of dynamic regulation of histone acetylation marks to control gene expression patterns associated with changes to EC behaviour.
We generated single-cell transcriptomic data using neonatal ECs to characterise the intra-tissue heterogeneity of the angiogenic retinal endothelium and studied the transcriptional changes resulting from deletion of Hbo1. Although we hypothesised that EC identity specification would be perturbed in the absence of HBO1, to our surprise we instead found that EC fate was mostly established normally. In our scRNA-seq dataset, we observed proliferating, arterial and venous/capillary ECs clusters that comprised both genotypes in equal proportions, suggesting no consistent change in their identity specification. This was in accordance with the normal morphological patterning of arteries and veins observed in retinas in which Hbo1 had been deleted. Interestingly, we found that tip cells were proportionally more abundant in retinas when Hbo1 was reduced. We were unable to distinguish stalk cells, possibly owing to overlapping gene expression patterns with arterial ECs, as both EC subtypes express Notch/TGFβ/BMP-activated genes, including Hey1, Id1 and Jag1 (Kalucka et al., 2020; Moya et al., 2012).
Our observation that tip cells were more abundant in the absence of HBO1 suggests that tip cells are the default EC type formed in the absence of HBO1-dependent chromatin organisation. This may be due to a role for HBO1 in normally restraining tip cell fate and/or regulating exit from the tip cell fate during angiogenesis. A role for HBO1 in restraining tip cell fate is supported by the finding that tip cell-expressed genes, including Esm1, Cxcr4 and Igfbp3, were upregulated when Hbo1 was deleted. We predict that this occurs through the actions of an intermediate transcriptional repressor that depends on H3K14ac for association at these gene loci. Both Esm1 and Igfbp3 are expressed by tip cells (del Toro et al., 2010; Rocha et al., 2014), but deletion of these factors does not cause a failure to make tip cells (Dallinga et al., 2020; Rocha et al., 2014). In contrast, deletion of Cxcr4 specifically from tip cells results in reduced likeliness of taking on the tip cell position (Pitulescu et al., 2017), suggesting that increased Cxcr4 expression in retinas in which Hbo1 has been deleted might contribute to the increase in tip cells. Alternatively, the defective migration observed in Hbo1 tip cells may result in prolonged exposure to the elevated VEGFA levels normally present at the sprout front and a subsequent increase in tip cell-associated genes. The increase in tip cells may also suggest a failure of Hbo1-deficient ECs to exit the tip cell phenotype. Classically, tip cell fate was thought to be laterally inhibited by activation of Notch signalling in stalk cells (Geudens and Gerhardt, 2011) or in a more inductive manner (Jakobsson et al., 2010; Bentley et al., 2014; Hasan et al., 2017; Pitulescu et al., 2017). Overactivation of Notch signalling leads to upregulation of the tip cell markers Esm1 and Cxcr4 (Pitulescu et al., 2017), but reduced sprouting and vessel growth (Pitulescu et al., 2017; Phng et al., 2009; Izumi et al., 2012), phenocopied by Hbo1-deficient retinas. Gene set testing did not show a significant enrichment for upregulation of Notch signalling genes.
The reduction in vascular outgrowth and vessel sprouts in retinas in which Hbo1 has been deleted may be caused by a defect in directed tip cell migration. Tip cells establish front-rear polarity during migration (Dubrac et al., 2016; Franco et al., 2015; Moya et al., 2012; Kim et al., 2019), thus we used front-rear polarity as a readout for an EC migration defect. In the absence of HBO1, tip cell front-rear polarity was disorganised relative to the direction of sprout elongation, suggesting a failure of biased migration that was consistent with the observed reduction in vessel growth. This was not due to defective FOXO1 signalling downstream of hypoxia, as shown by others (Kim et al., 2019). The Odc1 gene, which encodes ornithine decarboxylase, a key enzyme involved in polyamine biosynthesis (Pegg, 2006), was one of two genes downregulated in tip cells. Ornithine decarboxylase activity has been reported to promote in vitro vessel sprouting and EC migration through regulation of actin cytoskeleton dynamics (Kucharzewska et al., 2010), raising the possibility that its dysregulation is one cause of the Hbo1 migration phenotype. The observation that, although the number of vessel sprouts was reduced, the average sprout length was not different between control and Hbo1iEC/iEC mice suggests that tip cells that become front-polarised are able to migrate normally. Therefore, our data suggest that directionless polarity decreases the efficiency by which the ECs can migrate and expand the vessel network. Neovessel formation in the OIR model is dependent on directed tip cell migration (Dubrac et al., 2016). Therefore, our findings that neovascular lesions were reduced in the absence of HBO1 are also likely to be explained by inefficient tip cell migration.
The role of HBO1 in promoting sprouting EC behaviour is in accordance with the findings of a previous study that assessed HBO1 knockdown in HUVECs and zebrafish (Yan et al., 2018). However, our results do not support their conclusion that this was due to dysregulated VEGF signalling. We found that ECs in Hbo1iEC/iEC retinas expressed normal levels of Vegfr2. Additionally, Hbo1iEC/iEC mice were able to upregulate VEGFA-responsive ESM1 and DLL4, suggesting that at least some VEGFR2 signalling outputs are normal. Owing to the small number of tip cells in our dataset, we could not examine gene expression in the tip cell cluster with statistical robustness and therefore cannot rule out the possibility of altered VEGFA signalling pathway genes specifically in this cluster. Instead, based on the accumulation of ECs in the sprouting front and increased relative tip cell abundance when HBO1 is absent, we suspect that the number of tip cells present, together with a migration defect, dictates sprouting efficiency.
Although several studies have suggested a role for HBO1 in regulating proliferation and DNA replication licensing (Iizuka and Stillman, 1999; Johmura et al., 2008; Wu and Liu, 2008), this role is disputed. We have previously shown that proliferation in mouse embryonic fibroblasts is not affected by loss of HBO1 (Kueh et al., 2011) and human MCF7, HEK293T and HeLa cells can proliferate without HBO1 (Kueh et al., 2020). Consistent with other in vivo studies of HBO1 (Kueh et al., 2011), we did not find HBO1 to be essential for endothelial proliferation by EdU labelling. Supporting this, our scRNA-seq data did not reveal any differences in proliferating cells or regulation of genes relating to proliferation.
The relative abundance of upregulated genes when Hbo1 is deleted, as shown here in ECs and previously in HeLa cells (Kueh et al., 2020), suggests that the presence of H3K14ac at target genes is not simply for gene expression activation. Given the extensive occupancy of H3K14ac across the genome, it is possible that this mark is important for establishing widespread chromatin organisation to facilitate gene expression programmes. The requirement for HBO1 as an organiser of chromatin architecture is supported by time-course analysis of cell lines in which HBO1 has been deleted that show that loss of other histone acetylation marks occurs secondarily to loss of H3K14ac (Kueh et al., 2020). These observations suggest that loss of H3K14ac leads to a complex array of primary and secondary effects affecting multiple chromatin modifications across the chromatin landscape rather than HBO1 being directly responsible for the expression of a discrete set of target genes.
Our analysis of GO terms relating to genes differentially expressed in Hbo1iEC/iEC ECs provided some insights into the programmes that H3K14ac regulates. Strikingly, we found that many of the upregulated GO terms were associated with immature vessel phenotypes and behaviour, such as ‘blood vessel development’, ‘angiogenesis’ and ‘vasculature development’. Therefore, although EC identity was specified normally, they were maintained in a more immature state. This finding is consistent with our hypothesis that HBO1 would be required when new expression programmes are activated. We predict that the genes that are affected are those regulated by factors that are sensitive to loss of H3K14ac. Our finding that the adult vasculature did not express H3K14ac suggests that HBO1-dependent H3K14ac is unlikely to be required for maintaining gene expression patterns. Rather, this supports the proposed role for HBO1-dependent H3K14ac in organising the chromatin architecture required for new gene expression programmes to be established, such as during neonatal sprouting angiogenesis. The surprising finding that genes involved in pathways such as ‘angiogenesis’ and ‘blood vessel development’ were upregulated in the absence of Hbo1, despite a reduced angiogenesis phenotype is likely explained by a disturbance in the number of angiogenic tip cells and a migration defect. The outcome is more tip cells that highly express angiogenic genes but fail to migrate and grow the vasculature normally.
Overall, our data suggest that HBO1 is necessary for widespread chromatin organisation through its acetylation activity at H3K14 that enables ECs to activate new gene expression programmes. This is important in settings such as during vessel maturation and tip cell-directed normal and pathological vessel growth, but is redundant in settings of vascular quiescence. scRNA-seq highlighted the pathways that are dependent upon the activity of H3K14ac and suggested that tip cells are the default cell type formed when this transcriptional organisation is disrupted.
MATERIALS AND METHODS
Conditional Hbo1 mice (Kueh et al., 2011), Cdh5(PAC)-creERT2 mice (Wang et al., 2010) and ROSA26Sortm4(ACTB-tdTomato,-EGFP)Luo (mTmG) (Muzumdar et al., 2007) have all been previously described. Animals were maintained on an inbred C57BL/6 background. The day of birth was termed P0. Mice of both sexes were used. To induce Cdh5(PAC)-creERT2, mice were injected with 50 µg tamoxifen [MP Biomedicals, dissolved in sterile corn oil (Sigma-Aldrich, C8267) plus 5% ethanol] by intragastric injection at P1 and P2 using a 0.3 ml syringe with attached 29G needle (BD, 326103). Six-week-old adult mice were treated by oral gavage with 150 mg/kg tamoxifen prepared as above for three consecutive days and analysed at 10 weeks. Control mice were the following genotypes: Hbo1fl/+;Cdh5(PAC)-creERT2+ (heterozygote following tamoxifen administration), Hbo1fl/+ and Hbo1fl/fl.
All experiments involving animals were performed with procedures approved by The Walter & Eliza Hall Institute of Medical Research Animal Ethics Committee.
Nursing dams with P7 pups were housed in a Perspex chamber (BioSpherix) and exposed continuously to 74±1% oxygen in air maintained by a ProOx110 oxygen controller (BioSpherix). Mice were exposed to high oxygen from P7 to P12 before returning to normal room air (∼20% oxygen). Pups were fostered to BALB/c females following exposure to 3 days of high oxygen to prevent oxygen toxicity in dams. Neovascularisation was assessed at P17.
For whole-mount immunohistochemistry, eyes were fixed for 2 h in 4% paraformaldehyde at 4°C before retinas were dissected and blocked for 1 h at room temperature in Dulbecco's phosphate buffered saline (DPBS) with 1% Triton X-100 and 2% donkey serum. Retinas were stained with primary antibodies prepared in blocking solution overnight at 4°C, washed in DPBS containing 0.01% Triton X-100, then stained overnight with secondary antibodies prepared in blocking solution. Primary antibodies were: rat anti-PECAM1/CD31 (BD Pharmingen, 553370, 1/100), goat anti-PECAM1/CD31 (R&D Systems, AF3628, 1/200), rabbit anti-H3K14ac (Millipore, 07-353, 1/500), goat anti-collagen IV (Merck, AB769, 1/100), rabbit anti-cleaved (active) caspase 3 (CST, 9664, 1/50), rat anti-VE cadherin (BD Pharmingen, 555289, 1/500), rabbit anti-FLI1 (Abcam, ab15289, 1/200), human anti-ANG2 (4H10, 1/200) (Han et al., 2016), goat anti-ESM1 (R&D Systems, AF1999, 1/100), goat anti-DLL4 (R&D Systems, AF1389, 1/100), rabbit anti-EGFP (Invitrogen, A11122, 1/500) and rabbit anti-GM130 conjugated to A647 (Abcam, ab195303, 1/500). Secondary antibodies were: donkey anti-rabbit-Cy3 (Jackson ImmunoResearch, 711-165-152), donkey anti-rabbit-AF647 (Jackson ImmunoResearch, 711-605-152, donkey anti-rat-DL488 (Jackson ImmunoResearch, 712-485-153), donkey anti-rat-Cy3 (Jackson ImmunoResearch, 712-165-150), donkey anti-rat-AF647 (Jackson ImmunoResearch, 712-605-150), donkey anti-goat-DL405 (Jackson ImmunoResearch, 705-475-147) and donkey anti-goat-Cy3 (Jackson ImmunoResearch, 705-165-147). All secondary antibodies were used at 1/400. Retinas were incubated in 2 µg/ml Hoechst 33342 (Invitrogen, H3570) in 0.01% Triton X-100 DPBS for 2 h at room temperature and imaged the same day. Retinas were mounted with Prolong Diamond (Invitrogen, P36961). EdU was dissolved in DMSO to a concentration of 10 mg/ml. P6 pups were administered 100 µg EdU diluted in DPBS by intraperitoneal injection using a 0.3 ml syringe with attached 29G needle (BD, 326103). Cells were labelled with EdU for 2 h before euthanisation by cervical transection. Eyes were dissected and stained as per the procedure described above. Click-iT EdU Alexa Fluor 647 (Invitrogen, C10340) was performed as per manufacturer's instructions prior to mounting retina onto glass slides.
Imaging and image analysis
Retinas were imaged using a Leica TCS SP8 confocal microscope using 10×/0.4 NA, 20×/0.75 NA or 40×/1.30 NA objectives and Leica Application Suite software. All image analysis was performed in the Fiji distribution of ImageJ software (Schindelin et al., 2012). H3K14ac expression was analysed as the mean pixel intensity within endothelial nuclei, identified by Hoechst-33342 staining within PECAM1 signal using individual z-slices. Entire EC nuclei were traced manually using the freehand selection tool analysing >30 EC nuclei per image. H3K14ac levels in ECs were normalised to mean H3K14ac expression in >30 non-EC nuclei (determined as being PECAM1 negative). Vessel area was calculated using PECAM1 signal from maximum intensity projection images following application of a median filter (2 pixels) and ‘despeckle’ filter prior to manually adjusting the threshold and measuring the area. Radial outgrowth was quantified as mean distance from the optic nerve head at the centre of the retina to the edge of the vasculature. Branch points and segment density were counted within two regions per retina from maximum projection images of collagen IV staining (P6 retinas) and PECAM1 staining (adult retinas). P6 regions contained some arterial and venous plexus as well as some of the sprouting front plexus and adult regions were from near an artery in the peripheral retina. These images were despeckled, and a 2-pixel median filter was applied and made binary by thresholding area. The images were then skeletonised and the skeleton was analysed. Results output included details of the number of junctions (branch points) as well as the number of segments. These results were normalised to the area of the cropped region analysed. The proportion of retina vascularised by superficial and deep vessel layers was calculated by manually tracing around the vascularised retina and normalising to the manually traced total retina area. Apoptotic ECs [defined as cleaved (active) caspase 3+/PECAM1+ cells enclosed by collagen IV signal] were quantified manually from confocal z-stack images. Vessel regression was determined as a ratio of PECAM1+ vessel segment length to collagen IV+ vessel segment length. Binary masks of both the PECAM1 and the collagen IV channel were made by various morphological filters and thresholding the signal manually. The collagen IV+ PECAM1− vessel segment mask (i.e. regressing vessels) was generated by subtracting the PECAM1 mask from the collagen IV mask. The collagen IV+ PECAM1− mask and collagen IV mask were then skeletonised and length of vessels within each mask was measured. The collagen IV+ PECAM1−:collagen IV ratio was generated automatically based on vessel length. Data are represented as the PECAM1:collagen IV ratio. Vessel regression was measured in multiple regions across the retina, including areas of the high-density sprouting front and less-dense remodelling vessels. Proliferation was determined by labelling cells with EdU and quantified from maximum projection retina overview images with the optic nerve head area removed. Analysis was automated and binary masks of PECAM1, FLI1 and EdU signals were generated using automatic thresholding and various morphological filters. The output generated total EC number across retina (defined as PECAM1+ FLI1+) as well as total proliferating cells across the retina. Proliferating ECs were defined as PECAM1+ FLI1+ EdU+. ESM1+ ECs were quantified by drawing a line along the P6 sprouting front and manually counting ESM1+ FLI1+ PECAM1+ cells within 150 µm of this line. To distinguish ESM1+ ECs, individual z-slices were examined. ESM1+ cells were normalised to the total length of retina analysed. ESM1 expression was measured in individual ECs by segmenting the vasculature into individual cells using PECAM1 and FLI1 masks and then determining ESM1-positive and -negative ECs by manually thresholding the ESM1 signal. Following this, the mean level of ESM1 signal in ESM1-positive and -negative cells was measured. Data are presented as mean ESM1 signal in ESM1+ cells divided by mean ESM1 signal in ESM1− cells (background). ANG2+ ECs were quantified in a similar manner, instead using Hoechst 33342 as a nuclear marker to distinguish ECs. Vessel sprouts were defined as filopodia rich EC-body projections that extended from the edge of the vasculature and were counted manually along the sprouting front in one half of the retina and normalised to the total length of retina analysed. To measure tip cell front-rear polarity, Golgi apparatus were visualised using GM130 and the orientation was analysed in ESM1+ ECs along the sprouting front in one half of the retina. The following criteria were used to assign Golgi apparatus orientation: (1) Golgi apparatus were defined as front-polarised if part of the Golgi apparatus was between the sprouting edge and nucleus; (2) Golgi apparatus were defined as rear-polarised if part of the Golgi apparatus was at the opposite side of nucleus to sprouting edge; and (3) Golgi apparatus were defined as middle-polarised if the Golgi apparatus was neither front- nor rear-polarised or if Golgi apparatus was both front- and rear-polarised. Golgi apparatus orientation was represented as the percentage in each position. For EC density and vessel width analysis, the number of FLI1+ ECs in each region was manually counted and normalised to vessel area for EC density analysis or to total vessel length for vessel width analysis. Total length of vessels was calculated by skeletonising the filtered PECAM1 signal and summing all vessel segments together. In the sprouting front, this was measured in two regions at the edge of vessel plexus around the end of a vein and an artery. In the remodelling zone, this was measured in two regions between an artery and vein. Neovascular and vaso-obliterated area were quantified based on collagen IV and PECAM1 signals, respectively, which were segmented manually in Adobe Photoshop CC 2015 and the area calculated in Fiji. Neovascular area per retina was normalised to total retina area. Where necessary, a despeckle filter was applied to select channels in images displayed in figures for clarity. All data are shown as mean± s.e.m. Statistical analyses were performed for all quantitative data using Prism 9.0 (GraphPad).
Quantitative real-time PCR on lung endothelial cells
P6 lungs were removed, cut into smaller pieces and dissociated using 0.26 WU/ml Liberase TM (Roche, 05401119001), 10 µg/ml DNase 1 (Worthington Biochemical Corporation, LS002139) in DPBS at 37°C for 40 min with constant agitation. Samples were washed with 2% fetal calf serum in KDS-BSS salt solution and filtered through a 70 µm filter mesh cap (Falcon). To deplete red blood cells, samples were incubated with rat anti-TER119 antibody (BD Pharmingen, 553671) for 20 min at 4°C. Cell number in each sample was automatically counted using a Countess cell counter (Invitrogen), and Dynabeads (Invitrogen, 11035) were added at a ratio of one bead:one cell and incubated with cells for 20 min at 4°C. Red blood cells were depleted by applying the sample to a magnet and transferring the supernatant to a new tube. Depleted samples were stained using anti-PECAM1-APC (eBioscience, 17-0311, 1/100), anti-ICAM2-FITC (BD Pharmingen, 557444, 1/500) and anti-CD45-PerCPCy5.5 (BD Pharmingen, 552950, 1/200). ECs were sorted using an ARIA IIu sorter (BD Bioscience) twice. ECs were defined as PECAM1hiICAM2+CD45.2−. Dead cells were excluded based on DAPI uptake. The first sort was using the ‘yield’ parameter with the collected sample spun down and resorted using the ‘4-way purity’ parameter. ECs were pelleted and frozen. RNA was isolated using a Qiagen RNeasy Microkit (Qiagen, 74004) according to the manufacturer's instructions and quantified using an Agilent Tapestation 2200. cDNA was synthesised using SuperScript II reverse transcriptase (Invitrogen). Real-time PCR was performed in triplicate using SensiMix SYBR Hi-ROX (Bioline, QT605-05)-based detection on a LightCycler 480 (Roche). Hbo1 gene expression was normalised to expression of Hsp90.
scRNA-seq of retinal endothelial cells
scRNA-seq was performed on ECs isolated from P6 retinas of pairs of littermate control and Hbo1iEC/iEC mice. Control genotype mice were Hbo1fl/+Cdh5-creERT2+. Retinal ECs from four Hbo1iEC/IEC mice and four control mice were sorted into 384-well plates, sequenced and analysed, with one control mouse used across two plates, hence Fig. S3A and Fig. 4B show five control samples. For all other analyses, the ECs from the control mouse on two plates were pooled and analysed as one sample. For ease of isolation of ECs from the retina, the mTmG allele was intercrossed with control and Hbo1iEC/iEC mice. Both retinas per animal were dissected and dissociated using 0.26 WU/ml Liberase TM, 10 µg/ml DNase 1 in DPBS at 37°C for 30 min with constant agitation. Samples were washed with 2% fetal calf serum in KDS-BSS salt solution and filtered through a 50 µm filter mesh (Sefar, 03-50/31) into a 4 ml tube and pelleted. Neural retinal cells were depleted by incubating with rat anti-CD24 (BioLegend, 101803, clone M1/69, 1/100) and rat anti-CD73 (BioLegend, 127203, clone TY/11.8, 1/100) antibodies for 20 min at 4°C. Samples were then incubated with 25 µl of Dynabeads Biotin Binder beads (Invitrogen, 11047) for 20 min at 4°C and depleted by applying the sample to a magnet and transferring the supernatant to a new tube. Depleted samples were stained with anti-PECAM1-PECy7 (eBioscience, 25-0311-82, clone 390, 1/100) for 20 min on ice. Single ECs were flow-sorted into a chilled 384-well PCR plate (Greiner, 785290) containing 1.2 μl of primer/lysis mix [20 nM indexed polydT primer (custom-made, IDT), 1:6,000,000 dilution of ERCC RNA spike-in mix (Ambion, 4456740), 1 mM dNTPs (NEB, N0446S), 1.2 units SUPERaseIN RNase Inhibitor (Thermo Fisher, AM2696), DEPC water (Thermo Fisher, AM9920)] using a BD FACSAria III flow cytometer (BD Biosciences). Prepared RNA mixture plates were sealed and stored at −80°C. ECs were defined as PECAM1hiEGFP+. Dead cells were excluded based on DAPI uptake. Between 171 and 185 single cells per animal were sorted. Apart from one plate that contained only control genotype ECs, all other plates contained ECs from pairs of littermate control and Hbo1iEC/iEC retinas. In total, five plates were collected, and processed for sequencing and downstream analysis. Single-cell transcriptome libraries were generated using the CEL-Seq2 protocol (Hashimshony et al., 2016) with the following adaptations. Second strand synthesis was performed using NEBNext Second Strand Synthesis module (NEB, E6111S) in a final reaction volume of 8 µl and NucleoMag NGS Clean-up and Size select magnetic beads (Macherey-Nagel, 7449970.5) were used for all DNA purification and size selection steps. The libraries were sequenced on an Illumina NextSeq500 sequencer using an Illumina 75-cycle sequencing kit (sequencing run set-up: Read 1–14 bp, Index read–6 bp, Read 2–72 bp).
CEL-Seq2 scRNA-sequencing reads were mapped to the GRCm38.p6 mouse genome and ERCC spike-in sequences using the Subread aligner (v2.2.4) (Liao et al., 2013) and assigned to genes using scPipe (v1.10.0) with GENCODE vM18 primary assembly annotation. Gene counts were exported as a matrix by scPipe with UMI-aware counting. All subsequent analyses were performed in R (version 4.0.0) (R Core Team, 2020) with Bioconductor (version 3.11) (Huber et al., 2015) and were adapted from ‘Orchestrating Single Cell Analysis with R’ (Amezquita et al., 2020). Owing to differences in quality control (QC) metrics by plate, outlier thresholds were calculated per plate. Cells were removed from further analysis if they failed to achieve QC cutoffs as shown in Table S6 for total read counts, total genes detected, percentage of reads from spike-in ERCC RNA and percentage of reads from mitochondrial genes. Overall, 229 cells were excluded, leaving 1393 for downstream analysis. All genes detected in at least one sample were retained for downstream analysis (33,655 genes). Size factors were estimated using the normalisation by the deconvolution method (Lun et al., 2016a), as implemented in the scran package (v1.16.0), and used to compute log-transformed normalised expression values (logcounts). We estimated the per-gene ‘biological’ variation as the total variation minus an estimate of the ‘technical’ variation (estimated by fitting a plate-specific, mean-dependent trend to the variance of the spike-in transcripts). All genes with positive ‘biological’ variation that were not annotated as pseudogenes, ribosomal protein genes, genes on the sex chromosomes, or mitochondrial genes were declared as highly variable genes (HVGs). We then performed principal component analysis (PCA) of the logcounts of the HVGs and performed further dimensionality reduction for the purposes of visualisation using the UMAP method (McInnes et al., 2018). These principal components were used to build a shared nearest neighbour graph (Xu and Su, 2015), which was itself used as input to the Louvain algorithm, as implemented in the igraph package (v1.2.5) (https://igraph.org/), to identify clusters of similar cells. We next corrected for plate-specific differences using the mutual nearest neighbours (MNN) method (Haghverdi et al., 2018) as implemented in the batchelor package (v1.4.0). The MNN-corrected data were used in place of the principal components and data were re-clustered and visualised by UMAP. Clusters with highly similar expression patterns and assumed to underlie the same EC phenotype were merged in the same cluster. Next, we identified top-ranking marker genes for each cluster by performing pairwise differential gene expression analysis for each cluster against all other clusters using a Welch t-test. The results were ordered based on the overall log2 fold change, with the highest log2 fold change representing cluster-specific markers. We used previously published known markers of EC phenotypes as indicated in the text and other single-cell endothelial transcriptomic analyses (Vanlandewijck et al., 2018; Kalucka et al., 2020; Goveia et al., 2020) to annotate the clusters we identified based on their top-ranking marker gene expression. Annotation of proliferating clusters was performed using the prediction method implemented in the ‘cyclone’ function from the scran package (v1.16.0) (Scialdone et al., 2015). To perform differential expression analyses, we first created pseudobulk samples by summing the single-cell counts from the same sample and, optionally, cluster (Lun and Marioni, 2017; Crowell et al., 2020 preprint). These pseudobulk samples were used as input for differential expression analyses using the quasi-likelihood methods from the edgeR package (v3.30.3) (Chen et al., 2016; Lun et al., 2016b). Heatmaps were generated on the pseudobulk logCPM values using the pheatmap package (v1.0.12) with row normalisation. GO term and KEGG pathway analysis were performed using ‘goana’ and ‘kegga’ functions from the limma package (Ritchie et al., 2015). Comparison of differentially expressed genes to other pathways involved in sprouting angiogenesis was assessed by rotation gene set testing using the ‘fry’ function from edgeR and visualised by barcode plots. Significance was determined by two-sided P-value testing. Differential abundance analysis was performed on samples normalised for total number of cells using the edgeR package.
HUVECs (Lonza, C2519A) were grown in 15 cm2 dishes with EGM-2 media (Lonza, CC-3162) under hypoxic (3% O2) conditions. Each 15 cm2 dish was considered sufficient for two immunoprecipitation preparations (IP preps). HUVECs of 90-100% confluence were stimulated with 25 ng/ml VEGF165 (Peprotech, 100-20) for 2 h, then fixed with 1% methanol-free formaldehyde for 10 min at room temperature and quenched with 0.125 M glycine. Fixed cells were scraped with DPBS, transferred to a tube and pelleted at 400 g for 5 min at 4°C. Cells were lysed using 1× Buffer A (Cell Signaling Technology, 14282S) supplemented with sodium butyrate, protease inhibitors and DTT on ice for 15 min. Nuclei were pelleted at 1500 g for 5 min at 4°C and resuspended in 1× Buffer B (Cell Signaling Technology, 14282S) supplemented with sodium butyrate, protease inhibitors and DTT. Resuspended nuclei were passed through a 100 µm filter mesh and centrifuged at 1500 g for 5 min at 4°C. Nuclei were resuspended in 100 µl supplemented Buffer B per IP prep and transferred to a 1 ml tube. Chromatin was fragmented to between 150 bp and 900 bp using 0.15 µl micrococcal nuclease (Cell Signaling Technology, 10011) per IP prep for 20 min at 37°C. The digest was stopped by adding 10 µl of 0.5 M EDTA per IP prep and placing on ice. Nuclei were pelleted at 16,000 g for 1 min at 4°C. The nuclear pellet was resuspended in 150 µl ChIP lysis buffer (50 mM Tris-HCl pH 8, 5 mM EDTA, 0.5% SDS, supplemented with protease inhibitors and sodium butyrate) per IP and incubated on ice for 10 min. Lysates were sonicated to total power of 60 J to break the nuclear membrane. Lysates were cleared by centrifuging at 16,000 g for 10 min at 4°C and the supernatant was transferred to a new tube. Concentration was determined using A280 Nanodrop reading. For each H3K14ac and panH3 IP 40 µg chromatin was used and 80 µg chromatin was used per IgG IP. Samples were diluted 1/10 with ChIP dilution buffer (20 mM Tris-HCl pH 8, 2 mM EDTA, 150 mM NaCl, 0.01% SDS, 1% Triton X-100, supplemented with protease inhibitors and sodium butyrate) in DNA LoBind tubes (Eppendorf, 0030108035) and 2% of the volume removed as input sample. Each IP prep was incubated with 15 µl Magna ChIP A+G beads (Merck, 16-663) for 1 h at 4°C and applied to a magnet to remove non-specific binding. Supernatant was transferred to a new LoBind tube and incubated overnight at 4°C with rabbit anti-H3K14ac (Abcam, ab52946, 5 µg per IP), mouse anti-pan H3 (Abcam, ab10799, 10 µg per IP) or rabbit anti-IgG (Abcam, ab172730, 10 µg per IP) antibodies. Samples were incubated with 15 µl pre-washed Magna ChIP A+G beads for 1 h at 4°C. Beads were washed three times with TSEI buffer (20 mM Tris-HCl pH 8, 2 mM EDTA, 150 mM NaCl, 0.1% SDS, 1% Triton X-100), and once each with TSEII buffer (20 mM Tris-HCl pH 8, 2 mM EDTA, 500 mM NaCl, 0.1% SDS, 1% Triton X-100), TSEIII buffer (10 mM Tris-HCl pH 8, 1 mM EDTA, 250 mM LiCl, 0.1% SDS, 1% sodium deoxycholate, 1% IGEPAL CA-630) and TE buffer (10 mM Tris-HCl pH 8, 1 mM EDTA). Chromatin was eluted by incubating beads with 75 µl ChIP elution buffer (1% SDS, 0.1 M NaHCO3) for 30 min at 65°C with constant agitation, applying the sample tube to a magnet and transferring to a new LoBind tube. Samples (including input) were reverse crosslinked with 20 µg RNase A at 37°C for 30 min and then 40 µg Proteinase K at 65°C overnight. DNA was purified using a Qiagen MinElute PCR Purification Kit (Qiagen, 28004) following the manufacturer's instructions. DNA concentration was determined using an Agilent Tapestation 2200. For each of three H3K14ac, two panH3 and two input samples, 5 ng of DNA per were used for sequencing library preparation, generated using an Illumina TruSeq DNA library prep kit. ChIP libraries were sequenced on an Illumina NextSeq 500 sequencer. H3K14ac samples were sequenced with twice the depth as panH3 and input samples, by loading more starting material.
ChIP-seq reads were aligned to hg38 using Rsubread (Liao et al., 2019). The ChIP-Seq libraries were GC corrected using Benjamini's and Speed's method implemented by the deepTools program (Ramírez et al., 2016). The genomic features – gene bodies, pseudogene bodies, TSS−1 kb to TSS+1 kb regions, TSS−2 kb to TSS regions, transcription end site (TES) to TES+2 kb regions, and intergenic regions defined as ≥10 kb away from any gene excluding blacklisted regions – were considered in the analysis using hg38 NCBI annotation. The number of read pairs overlapping each genomic region was summarised using featureCounts. Student's t-test was used to compare log2 counts per million (CPM) between the non-overlapping genomic regions. Differential expression analysis between H3K14ac and panH3 with respect to each of the genomic regions were undertaken using the edgeR and limma (Ritchie et al., 2015) software packages. Library sizes were normalised using the trimmed mean of M-values (TMM) method (Robinson and Oshlack, 2010). Differential expression was assessed using glmQLFit (quasi-likelihood generalised linear models) in edgeR (Chen et al., 2016). Over-representation of GO terms for the differentially expressed genes was identified using limma's goana function. Barcode plots illustrating the enrichment of gene sets in H3K14ac compared with panH3 were drawn using limma's barcodeplot function. Gene set enrichment tests used the roast method (Wu et al., 2010). The coverages for non-overlapping 5 bp bins for each ChIP-seq library were computed using the deepTools program.
The authors thank Ralf H. Adams (Max Planck Institute, Münster) for Cdh5(PAC)-creERT2 mice, Goh Young Koh (KAIST, Daejeon) for generously providing the anti-ANG2 antibody 4H10, Glenna-Faye Dabrowski, Lauren Wilkens, Andrea Briffa, Jaclyn Gilbert and Dannielle Cooper for expert animal care, and members of the Epigenetics and Development Division (Walter & Eliza Hall Institute) for informative discussion. This work was made possible through Victorian State Government Operational Infrastructure Support and Australian Government NHMRC IRIISS.
Conceptualization: Z.L.G., A.K.V., L.C.; Methodology: R.C.A.S., S.H.N.; Software: L.W.; Formal analysis: Z.L.G., P.F.H., W.A., A.L.G., A.K.V., L.C.; Investigation: Z.L.G., S.M.L., T.M.B.; Resources: R.C.A.S., T.T.; Writing - original draft: Z.L.G., A.K.V., L.C.; Writing - review & editing: Z.L.G., A.K.V., L.C.; Supervision: D.A.-Z., G.K.S., A.K.V., L.C.; Funding acquisition: A.K.V., L.C.
This work was supported by the National Health and Medical Research Council, Australia (Project Grants 1084248 and 1143612 to A.K.V. and T.T., and 1124812 to S.H.N.; Research Fellowships 1081421 to A.K.V., 1154970 to G.K.S. and 1003435 to T.T.; Investigator Grant 1176789 to A.K.V.); the Australian Research Council (Future Fellowships: 110100891 to L.C.); the L.E.W. Carty Charitable Fund (L.C.); an Australian Government Research Training Program Scholarship (to Z.L.G.); and the Royal Melbourne Hospital Neuroscience Foundation (R.C.A.S.).
The authors declare no competing or financial interests.