2D hPSC-ECs lose their endothelial identity during longer-term culture

hPSC-ECs were initially derived via a published 2D-adherent growth factor-based protocol [11]; this resulted in a primitive vascular network supported by non-ECs (Fig. 1A). CD31pos/NRP1pos hPSC-ECs were purified from the differentiation culture via fluorescence-activated cell sorting (FACS) and cultured in 2D, resulting in the adoption of the traditional endothelial cobblestone morphology (Fig. 1B), formation of tube-like structures in vitro (Fig. 1C), and expression of pan-endothelial marker proteins (Fig. 1D). The expression of endothelial- and mesenchymal-associated genes alongside endothelial transcription factors that are active during embryonic development was assessed by RT-qPCR at three time points during differentiation: day 5 (early differentiation), day 12 (pre-FACS), and day 19 (post-FACS expansion). Although CD31 expression was upregulated following FACS (Fig. S1A), the expression of mesenchymal markers remained consistently high (Fig. S1B), suggesting that the emerging 2D hPSC-ECs, and not just the supporting non-ECs, expressed mesenchymal-associated genes. Expression of ETV2, the master regulator of EC development, displayed a trend towards downregulation, whereas ERG and FLI1, both involved in endothelial lineage maintenance [17], and SoxF transcription factors were upregulated as the 2D differentiation ensued (Fig. S1C-D).

Despite the previously reported long-term in vitro stability of this CD31pos/NRP1pos double-positive subpopulation [11], elongated, mesenchymal-like cells began to emerge at passage 3 (Fig. 1Ei), disturbing the cobblestone morphology observed in freshly sorted cells (Fig. 1B). These mesenchymal-like cells proceeded to become the predominant cell type by passage 5 (Fig. 1Eii). This was driven by CD31pos/FSP1pos double-positive 2D hPSC-ECs (Fig. 1F) that were absent from both the starting culture and from passage-matched native ECs (Fig. 1G). It is plausible that these cells behave as an intermediary population of 2D hPSC-ECs undergoing EndMT. This was substantiated by the continual downregulation of both CD31 (Fig. S1E) and the pro-maintenance TGFβ receptor, ALK1 (Fig. S1F) in passaged 2D hPSC-ECs. Expression of TGFB2 remained largely unchanged (Fig. S1G); however, downregulation of ALKI may suggest a shift towards TGFβ2-mediated EndMT.

SOX7 is upregulated in 2D hPSC-ECs relative to native ECs

mRNA microarray analysis was conducted on 2D hPSC-ECs to ascertain their global gene expression and determine how this differed from native ECs (Fig. 2A). Owing to their longstanding use in endothelial research, human umbilical vein ECs (HUVECs) were utilised as the native population. CellNet analysis of the microarray data highlighted the transition from gene regulatory networks (GRNs) associated with hPSCs to those associated with ECs (Fig. 2B) following endothelial differentiation. GRNs for other mesodermal-derived cells (cardiac and musculoskeletal) were not activated. The mean classification score (MCS, the probability of the cells expressing genes at a level that is indistinguishable from their reference cell type) for hPSCs decreased during endothelial differentiation (Fig. 2C) as did gene expression of the pluripotency marker gene POU5F1 (Fig. 2D). Although activated, 2D hPSC-ECs had a lower MCS than HUVECs (0.72 versus 0.87, respectively. Figure 2E) suggesting only partial activation of the endothelial GRN. CellNet network analysis also identified candidate transcription factors differentially expressed between the four groups of cells (Fig. 2F). Subsequent RT-qPCR analysis substantiated the upregulation of seven genes, SOX7, SOX17, SOX18, YAP1, LYL1, HOXB7, and HOXB3 (herein, CellNet-identified genes), in D19 2D hPSC-ECs relative to HUVECs (Fig. 2G). With a view to fully recapitulating the native EC GRN, 2D hPSC-ECs were subjected to small interfering RNA (siRNA)-mediated downregulation of SOX7. This resulted in an augmentation of the CD31pos population (Fig. 2H), which was independent of an increase in the proliferation of CD31pos hPSC-ECs (Fig. 2I). These data suggest that SOX7 knock-out-dependent improvement in cellular composition was driven by enhanced cell stability and retention of the initial endothelial identity rather than instigating proliferation of the CD31pos population. Functional analysis, however, revealed that SOX7 siRNA-treated 2D hPSC-ECs generated less complex in vitro vascular networks comprised of unconnected tubes (Fig. 2J, K). Targeting SOX7 during the 2D endothelial differentiation protocol (Fig. 2L) did not affect the relative abundance of 2D hPSC-ECs (Fig. 2M) but rather impinged upon the differentiation process, resulting in downregulation of both CD31 and VE-cadherin (Fig. 2N) thereby suggesting less mature ECs that would therefore not be suitable for either disease modelling or CTE applications.

Matrigel-free, 3D suspension culture generates vascular organoids at scale

To account for the effects of the 3D microenvironment on the emerging hPSC-ECs, hPSCs were grown as spheroids in an adapted 3D suspension culture protocol devoid of Matrigel. CD31pos cells emerged within these structures by day 5 of the differentiation protocol, forming early vascular networks in contrast to age-matched 2D differentiation (Fig. 3A, B) and continued to develop in complexity throughout the duration of the differentiation protocol. As the CD31pos/NRP1pos subpopulation of 2D hPSC-ECs did not confer greater phenotypic stability (Fig. 1F), the total population of CD31pos cells was isolated via FACS and analysed (Fig. 3C, D). Endothelial differentiation of hPSCs in 3D suspension augmented the total number of CD31pos cells with a trend towards a greater yield (P = 0.06, Fig. 3E) and generated cells with greater size homogeneity, as assessed via forward scatter analysis (Fig. 3F, G), potentially indicating a more homogenous population of ECs. Despite the absence of Matrigel, 2D- and 3D-hPSC-ECs had comparable gene expression of the pan-endothelial markers (Fig. 3H, I); however, VEGFA expression was downregulated in 3D hPSC-ECs (Fig. 3J).

High concentration VEGF-A treatment facilitates long-term culture stability

The effects of VEGF-A supplementation on 3D hPSC-EC phenotype and function were evaluated. hPSC-ECs emerging from the 3D differentiation were cultured in either normal EGM2 (herein, 3D hPSC-ECs) or EGM2 supplemented with variable (0, 5, 10, 50 ng/ml) concentrations of VEGF-A (Fig. 4A). Cultures exposed to two weeks of normal EGM2 contained an abundance of mesenchymal cells; however, cultures subjected to increasing VEGF-A concentrations displayed a positive concentration-response to stability, with the highest dose (50 ng/ml) generating homogenous cultures of CD31pos 3D hPSC-ECs (herein, 3DV hPSC-ECs; Fig. 4B). High content image analysis demonstrated a 2.5-fold increase in the mean CD31pos area following high VEGF-A treatment; however, this did not translate to a greater degree of in vitro tubular structures which remained unchanged (Fig. 4C). RT-qPCR analysis revealed 3DV hPSC-ECs had an augmented expression of pan-endothelial markers in addition to SOX7, with endogenous VEGF-A displaying a trend towards upregulation (Fig. 4D). Moreover, administering 50 ng/ml VEGF-A to 3D hPSC-ECs that had undergone EndMT (Fig. S2A) resulted in the re-emergence of CD31pos cells (Fig. S2B) in these transdifferentiated cultures alongside the restoration of CD31 and PNP gene expression to levels comparable with freshly isolated 3D hPSC-ECs (Fig. S2C-E).

Fig. 1

hPSC-ECs attained via a 2D endothelial differentiation protocol results in a heterogeneous population: (a) 3D rendering of day 12 2D hPSC-ECs (CD31, green) generating a primitive vascular network supported by non-endothelial cell populations (FSP1, red). Scale bar represents 100 μm. Phase-contrast microscopy images of day 19 2D hPSC-ECs demonstrating (b) the classical endothelial cobblestone morphology in adherent monolayer culture and (c) in vitro tube formation capacity on growth factor-reduced Matrigel. Scale bar represents 500 μm. (d) Representative immunofluorescent images of day 19 2D hPSC-ECs stained for (i) CD31, (ii) VE-Cadherin, and (iii) von Willebrand factor, vWF. Scale bar represents 10 μm. (e) Phase-contrast microscopy images of (i) passage 3 and (ii) passage 5 2D hPSC-ECs displaying the emergence of mesenchymal-like cells. Representative immunofluorescent images of P3 (f) 2D hPSC-ECs and (G) HUVECs cultured for 7 days. Cells were stained for CD31 (green), FSP1 (red), and Hoechst (blue). The white arrows highlight the emergence of the 2D hPSC-EC CD31pos/FSP1pos cells undergoing EndMT. Scale bars represent 50 μm

3DV hPSC-ECs demonstrate an angiogenic protein profile comparable to that of primary cardiac HMVECs

Comparison of the different hPSC-ECs (2D-, 3D-, 3DV-hPSC-ECs) populations was first conducted at the protein level using a Human Angiogenesis Array immunoassay to determine the expression of 54 angiogenesis-associated proteins. As ECs derived from different vascular networks have unique protein expression profiles, two native EC populations, cardiac HMVECs (HMVEC-Cs) and human coronary artery endothelial cells (HCAECs), representative of the micro- and macro-vasculature, respectively, were used as references of native ECs (Fig. 4E).

Despite receiving high-concentration VEGF-A treatment, the 3DV hPSC-ECs exhibited the highest protein levels of endogenous VEGF-A and angiopoietin-2 (Fig. 4E). This differed from the 3D hPSC-ECs that instead expressed the angiogenesis inhibitor, thrombospondin 1, similarly expressed in HCAECs. Conversely, expression of tissue factor, the key activator of coagulation also implicated in angiogenesis [18], was limited to the groups derived under 3D conditions, indicating that its protein expression was induced independent of exogenous VEGF-A and rather due to signals from the 3D microenvironment. Expression of the endothelial dysfunction-associated protein, Pentraxin 3 [19], was restricted to the native EC populations. Endothelin-1, also implicated in endothelial dysfunction [20], was detected in all EC populations (highest in HCAECs) apart from the 3DV hPSC-ECs.

String protein network analysis of the 3DV hPSC-ECs identified activation of VEGF-related networks with VEGF-A acting as the central mediator (Fig. 4F). Overall, principal component analysis (PCA) conducted for the angiogenic-protein profiles of all EC populations revealed that the 3DV hPSC-ECs exhibited an expression profile similar to that of native HMVEC-Cs (Fig. 4G).

3DV hPSC-ECs express markers associated with cardiac microvascular ECs

Further in-depth characterisation of the different endothelial populations was conducted via scRNA-seq. Despite an apparent similarity in their angiogenic profile, unsupervised clustering following scRNA-seq unveiled that the 3D-derived cells formed two clusters independent from the HMVEC-Cs that were demarcated by their expression of CD31 (Fig. 5A, B). The CD31pos cells were considered hPSC-ECs, whilst the CD31neg cells potentially represented the mural cells contained within the vascular organoids. The 2D hPSC-ECs also clustered separately from both the HMVEC-Cs and the 3D-derived cells.

The gene expression of ECs arising from different subtypes, vascular beds, and organs has recently been characterised via scRNA-seq studies [21,22,23,24,25], resulting in the identification of marker genes indicative of diverse endothelial parameters that were evaluated herein. The robust expression of pan-endothelial markers by the 3DV hPSC-ECs substantiated their highly endothelial nature, a stark contrast to the comparatively low expression observed in the HMVEC-Cs (Fig. 5C). Moreover, the 3DV hPSC-ECs preferentially expressed vascular endothelial markers whilst robust expression of lymphatic [25, 26] and valvular [14] endothelial marker genes was identified in the 2D hPSC-ECs and HMVEC-Cs, respectively (Fig. 5C). In addition to SOX7 expression (Fig. 5C), the 3DV hPSC-ECs expressed several of the CellNet-identified genes at levels greater than the 2D hPSC-ECs (Fig. S3).

In line with their valvular profile (Fig. 5C) and expression of endothelial dysfunction-associated proteins (Fig. 4E), HMVEC-Cs also exhibited prominent expression of genes implicated in EndMT [27,28,29,30] (Fig. 5D). Notably, the EndMT-inducing transcription factors SNAI1 and BMP6 were also expressed by 2D hPSC-ECs whereas 3DV hPSC-ECs expressed genes associated with the inhibition of EndMT [31,32,33,34] (Fig. 5D). Physiological function was also assessed by investigating VEGF-associated genes (Fig. 5D). The VEGF-A protein expression observed in the 3DV hPSC-ECs (Fig. 4E) was corroborated by VEGFA gene expression. The 2D hPSC-ECs and the HMVEC-Cs, however, expressed VEGFB, associated with pathological angiogenesis [35], with the latter also expressing the lymphangiogenic VEGFC and VEGFD [36]. Although lacking VEGFC expression, 3DV hPSC-ECs expressed its receptor, FLT4, in addition to the VEGF-A receptors, FLT1 and KDR.

The vascular bed(s) to which the different EC samples belonged to was then evaluated. As projections of the arteries, arterioles express artery-associated genes in addition to their own specific markers [22, 25, 37] (Fig. 5E). 3DV hPSC-ECs had the highest expression of this gene set with the notable exception of MECOM, which is recognised for its role in capillary arterialisation [21], and was instead expressed highly by the HMVEC-Cs. In contrast, SMAD1 is required for the formation of capillaries [21] and was highly expressed by the 3DV hPSC-ECs alongside capillary-associated genes [24, 25, 38, 39] (Fig. 5E). Conversely, HMVEC-Cs were devoid of these arteriole and capillary markers, thus further questioning their microvascular identity.

Investigation of venous markers revealed preferential expression by the 2D hPSC-ECs and HMVEC-Cs (Fig. 5E). Lymphatic endothelial cells (LECs) are derived from venous ECs during development [40] and therefore share canonical venous markers, including NR2F2, EPHB4, and NRP2 [25, 41, 42]. In contrast to arterioles, fewer venule-specific markers have been established. While the chemokine receptor ACKR1 has been reported in murine venules [43], its expression was largely absent from all EC populations (Fig. 5E). High endothelial venules (HEV) are specialised post-capillary venules that facilitate the migration of lymphocytes from the bloodstream to lymphatic vessels [44]. Consistent with their lymphatic association (Fig. 5C), 2D hPSC-ECs and HMVEC-Cs also expressed genes associated with the HEVs [45] (Fig. 5E).

In accordance with their capillary-like gene expression profile (Fig. 5E), the 3DV hPSC-ECs were highly enriched for genes associated with tip cells (Fig. 5E) that guide the angiogenic sprout in the direction of the angiogenic stimuli. Tip cells are identified through several markers, including CD34 [46], and encode genes associated with augmenting angiogenesis [47], metalloproteases for degrading the basement membrane [47, 48], collagen for subsequent basement membrane reconstruction [49], and genes responsible for resolving angiogenic sprouts to prevent aberrant angiogenesis [50, 51]. HMVEC-Cs did not express this gene set, and whilst a subset of these genes exhibited low expression in the 2D- and 3D-hPSC-ECs, the complete profile was robustly detected in the 3DV hPSC-ECs.

In a cardiac context, cardiac microvasculature is essential beyond angiogenesis to maintain cardiac function. Therefore, genes identified in cardiac ECs that mediate such functions were also investigated [52,53,54,55,56,57,58]. Indeed, the 3DV hPSC-ECs not only expressed genes associated with coronary angiogenesis [52, 53] but also with cardiomyocyte metabolism [54] and the prevention of EndMT [55, 56]. Further, cardiac ECs also exhibited modest expression of cardiomyocyte myofibrillar genes (CMFs) [59] that were heterogeneously expressed by both 3D-derived hPSC-ECs.

Although the 3DV hPSC-ECs expressed microvascular-associated genes (Fig. 5E) and genes that have been reported to facilitate cardiac function (Fig. 5F), further comparison with bona fide human cardiac ECs was required to confirm that high-concentration VEGF-A treatment mediated specification of the 3D hPSC-ECs towards a cardiac microvascular-like profile. The 3DV hPSC-ECs demonstrated robust expression of genes associated with cardiac capillary ECs (Fig. S4A) and cardiac arterioles (Fig. S4B), as identified in healthy adult human cardiac samples [60] and human foetal cardiac ECs [21, 61]. This was in contrast to genes associated with the cardiac macrovessels (Fig. S4C, D) and endocardial cells (Fig. S4E). Transcription factors that are preferentially expressed by different cardiac endothelial subtypes have also been identified from these human samples. To this end, the 3DV hPSC-ECs displayed a profile consistent with cardiac capillary- and microvascular-identified transcription factors (Fig. S4F) yet lacked expression of the cardiac arteriovenous transcription factors (Fig. S4G).

Upon integration with the human samples, the CD31pos component of the 3DV-derived cells clustered alongside the CD31pos endothelial cells emerging from the human samples whilst displaying divergence from the LYVE1 expressing cells that were indicative of the human cardiac lymphatic ECs (Fig. S5A-D). The CD31neg component of the 3DV-derived cells once again clustered separately (Fig. 5A, Fig. S5A, C) and exhibited pronounced expression of the conventional pericyte marker, PDGFRB (Fig. S5E) but not the fibroblast-specific marker, LUM (Fig. S6F) that was predominantly derived from the adult human samples. This further strengthened the notion that vascular organoids also generated mural cells that are a critical component of the microvasculature.

Owing to the abundance of non-endothelial cells within this integrated data set (Fig. S5C-E), the ECs were selected for further analysis through their expression of CD31 and the published data sets were downsampled to attain cell numbers comparable for all samples. The resultant CD31high cells from the 3DV hPSC-EC sample (Fig.