The human breast cancer cell lines were maintained in standard growth media (ATCC: MCF-7, T47D, SKBR3, MDA231) in RPMI (Gibco) plus 10 % FBS (Sigma) and PenStrep (Gibco); Asterland SUM149, SUM159 in Ham’s F-12 (Gibco) plus 5 % FBS (Sigma), insulin (Gibco 5ug/mL), and hydrocortisone (Gibco 1ug/mL). Human endothelial cells (EC) were all purchased from Lonza and were grown in EBM-2 media with BulletKit additives (Lonza CC-3202). Blood microvascular endothelial cells (BEC) and lymphatic microvascular endothelial cells (LEC) were from the same donor, and along with HUVECs, were all used within 5 passages. All vascular cell lines were grown in EBM-2 for 48 h prior to RNA extraction (monocultures and cocultures). Three-dimensional morphology experiments were performed by coating a Lab-Tek 8-well chamber slide (Thermo Scientific) with 125 uL of Matrigel (Becton–Dickinson) and then plating 50,000 cells in each well. To allow for discrimination of cancer cells and ECs, they were labeled with Sigma’s PKH67 (green) and PKH26 (red) dyes, respectively, prior to co-culture. All morphological studies (Fig. 4C–E) were performed for 18 h. Images of cell culture experiments were taken with a Nikon inverted phase contrast microscope and recorded with OpenLab software (Fig. 4; Supplemental Fig. 2). Confocal images (Fig. 4E) were taken with an Olympus FV 500 Confocal Laser Scanning Microscope and processed with Olympus FluoView software. Immunofluorescence images of xenografts tumors were acquired with an Olympus IX81 Inverted Light Microscope. All fluorescence images were combined with Image J v1.46.

RNA was prepared from human breast cancer cell lines and ECs with Qiagen’s RNeasy mini kit. Gene expression microarrays were performed according to established protocols [10, 11], with all microarray data publicly available at the UNC microarray database (UNCMD) https://genome.unc.edu/. New microarrays have been deposited in the Gene Expression Omnibus under the accession number GSE37145, with previously published data available under GSE31870. Prior to analyses, the expression data were downloaded from the UNCMD, and the probes were filtered by requiring the Lowess normalized intensity values in both sample and control to be greater than 10 dpi and present on more than 70 % of microarrays. The normalized log₂ ratios (Cy5 sample/Cy3 control) of probes mapping to the same Entrez gene ID were averaged and median centered to obtain the final dataset.

For the vascular content signature, RNA was prepared from human breast cancer cell lines (MCF7, T47D, SKBR3, SUM149, SUM159, MDA231), endothelial cells (HUVEC, LEC, BEC) and commercial RNA for human tissue was obtained from Ambion (brain; AM6050), (lung; AM7968), (liver; AM7960), (lymph node; AM7894), Clontech (bone marrow; 636591), and Biochain (bone marrow; R1234024-10). The breast cancer cell line arrays and the human organ arrays were compared against endothelial cell arrays from BEC, LEC, and HUVEC in a two-class SAM. Genes with high expression were used to generate this signature.

The activated endothelium signature was generated from microarrays of BEC and LEC RNA that had been grown in monoculture or transwell co-cultured with cancer cells. For transwell co-culture, 200,000 endothelial cells were plated into a well of a 6-well plate, then a transwell filter with 0.4 micron pores (Corning) was inserted into the wells and 200,000 cancer cells were added to the upper compartment in EBM-2. After 48 h the transwell inserts were removed and the endothelial cell RNA was extracted (Qiagen). For the activated endothelium signature a two-class SAM was performed on triplicate BEC and LEC arrays compared to arrays from BECs and LECs that had each been transwell cultured with the six breast cancer cell lines described above. Each tumor’s endothelial signature score was determined by averaging the log₂ expression values for all genes in the signature (either 74 for vascular content, or 110 for activated endothelium) that were also found in the different test datasets. To separate out the proliferation component of the signature, all genes with a Pearson correlation value greater than 0.5 to a 11-gene proliferation signature [12] were considered proliferation related; due to the reduced number of genes present, this signature was not able to be divided in the merged 550 tumor dataset. Lastly, a mouse mammary tumor gene expression dataset was also examined that has been previously published [13] (GSE3165 and GSE27101).

All statistical tests were performed with WinSTAT, R v2.15.1, and Cluster v3.0.

All animal procedures were done under a protocol approved by the University of North Carolina Animal Care and Use Committee. To establish MCF-7, MDA-231, SUM159 tumors, 3 Nod scid gamma (NOD.Cg-Prkdc^scid IL2rg^tm1Wjl/SzJ JAX^®) mice for each tumor type were anesthetized with isoflurane, and one-million cells in 100 % Matrigel with growth factors (Becton–Dickinson) were injected into the lactiferous duct of the fourth (inguinal) mammary gland. For MCF-7 tumors, mice were also implanted with an estradiol-releasing silastic pellet as previously described [14]. MDA-231 and SUM159 tumors were grown for 18 days, MCF-7 tumors were grown for 24 days. The difference in growth times reflects the amount of time needed to extract similar sized tumors. All tumors were removed when they were ~7 × 7 mm. To label perfusing vasculature in vivo, mice were injected with 1 mg of Texas red-conjugated dextran (molecular weight 70,000; Invitrogen/Molecular Probes, Eugene, OR) diluted in PBS (5 mg/mL) and then euthanized 5 min later, as previously described [15]. Prior to being embedded and frozen in Optimal Cutting Temperature (OCT) Compound, tissues were fixed in 10 % formalin overnight, and then 30 % sucrose overnight. Tissues were then stored at −80 °C until they were cut into 9–10 micron thick sections. Primary antibodies utilized were from Dako (vWF (A0082)), Novus Biosciences (LYVE1 [NBP1-43411], vimentin [NBP2-12472], cytokeratin 19 [NB100-79916]), Novacastra (CD34 [NCL-L-END]), and Santa Cruz (PECAM [sc-101454]). Secondary antibodies were from Molecular Probes/Life Technologies (Goat anti-rabbit 488) and Jackson ImmunoResearch (donkey anti-rat FITC). Mounting media containing DAPI was from Invitrogen (P36931).

The majority of nearly one-hundred publications have found that high microvessel density is associated with poor prognosis in breast cancer [16]. Since the results from this assay can vary depending on the vascular antibody utilized and/or one’s definition of a vascular hot-spot, we aimed to determine if vascular gene expression signatures could serve as an alternative biomarker to identify patients with an increased likelihood of distant metastasis or death. To test this hypothesis, we developed two novel Endothelial Cell (EC)-derived signatures and also tested a published EC signature that has been shown to be specific for the microvasculature [17]. For comparisons of these signatures, we used five previously published microarray data sets of breast cancer patients; these analyses comprised >3,000 human breast tumors with ~10 % overlap [2, 9, 18–20]. To contrast how the vascular signatures were expressed in normal breast samples, breast tumors, and pure endothelial cell lines, we determined average gene expression signature scores for each sample and in three EC lines (HUVEC, BEC, and LEC). To visualize these values across different subtypes of human breast tumors, normal breast samples, and EC lines, each sample’s signature was plotted as box-and-whisker plots (Fig. 1). On average, expression of the Wallgard et al. vascular signature was highest in EC lines, followed by normal breast reduction mammoplasty tissues (Fig. 1A). Interestingly, the claudin-low tumors had the highest vasculature signature expression when compared against any of the other tumor subtypes (t test, p < 0.0001). When this signature was examined on a database of mouse mammary tumors and normal mammary glands [13], the normal mammary tissue and murine claudin-low tumors also exhibited high expression (Fig. 1A).Fig. 1

Vascular gene expression signatures in different intrinsic subtypes of human breast tumors and transgenic mouse models of mammary cancer. Box-and-whisker plots are shown for five human breast tumor datasets (Left) and a mouse mammary tumor dataset [13] (Right). Gene expression signature scores were identified for endothelial cell lines (BEC, LEC, HUVEC), each breast cancer dataset, and then combined for display in the following order within each subtype; Combined 855 [9], MDACC [20], Merged 550 [19], METABRIC [18], UNC [2]. a Wallgard et al., b vascular content, c activated endothelium, and d VEGF/hypoxia signatures. The log₂ mean signature expression for each tumor is shown as a cross, the bar indicates the median value, whiskers show the range within subtype and are the 1.5 * inter-quartile range

We next sought to directly compare each endothelial signature to several other vascular signatures that have been previously identified [17, 23–25]. Pearson correlations were determined between the various signatures and a known cell proliferation signature [12] using the same five breast cancer datasets. Other than the Tumor Vascular A signature [24], as a whole, most vascular signatures were positively correlated (Fig. 2). We also found that the activated endothelium signature showed a strong correlation with proliferation (0.57–0.71). We therefore ‘separated’ the proliferation component from the rest of the activated endothelium signature by identifying genes with a Pearson correlation greater than 0.5 to the proliferation signature [26]. This resulted in two distinct signatures: an ‘activated endothelium proliferation component’ and an ‘activated endothelium non-proliferation component’ (Supplemental Fig. 1). In each of the datasets tested, the vascular content, Wallgard et al., and vasculogenic mimicry signatures were all correlated (>0.5), suggesting that these three signatures were tracking similar biological processes. The VEGF/hypoxia signature showed positive correlations with activated endothelium; interestingly, both the proliferative and non-proliferative components of the activated endothelium signature had smaller Pearson correlations with the VEGF/hypoxia signature than the complete activated endothelium signature, indicating that the VEGF/hypoxia signature identifies both of these biological processes.Fig. 2

Assessment of the relatedness of vascular gene expression signatures. Shown are Pearson correlation coefficients of gene expression signatures from the five breast cancer datasets. Positive values are colored red and negative values are green

To determine if the expression of any of the vascular signatures correlated with increased metastasis in vivo, univariate Cox proportional hazards models analyses were performed for each signature on the five datasets (Supplemental Table 3). Several variables individually were able to significantly predict metastatic relapse in every dataset, including the luminal B and basal-like subtype status (as compared to luminal A), ER status, activated endothelium (including both the proliferation and non-proliferation components), VEGF/hypoxia, vasculogenic mimicry, and the 11-gene proliferation signature. We next used multivariate Cox proportional hazards models to determine if any of the vascular signatures provided additional prognostic information in addition to intrinsic subtype classification. In these analyses, each vascular signature was individually tested against the metastasis predicting information contained within intrinsic subtype status; note that only the vascular signature information is shown in Table 1. Interestingly, both the activated endothelium and VEGF/hypoxia signatures were the only two signatures that were found to be significant in all datasets. Therefore these signatures can be used to identify particularly aggressive subsets of tumors within a given intrinsic subtype (the full table including intrinsic subtype information is shown in Supplemental Table 4). These findings support the conclusion that the quality of endothelium and the specific heterotypic interactions between endothelium and epithelium are both important components of metastatic progression.Table 1

Multivariate analyses of the vascular signatures tested individually against intrinsic subtype for metastasis prediction; tested in the five breast tumor datasets. Bolded signatures added metastasis predicting information to intrinsic subtype classification. Full table is presented in Supplemental Table 4

We next aimed to understand if any of the vascular gene expression signatures were a surrogate for the histology-based microvessel density (MVD) assay. Therefore, MVD scores were determined on a 70-tumor dataset that had also been subjected to gene expression microarrays (UNC70) (Supplemental Table 5). Although the sample set was small, on average MVD was similar across the breast cancer subtypes, except in the normal-like tumors which are comprised mostly of normal breast tissues (Fig. 3A). In this dataset, high MVD was significantly associated with decreased relapse free survival in a univariate analysis when tested as a continuous variable (p = 0.04), and was also trending towards significance in Kaplan–Meier analyses when the sample set was divided into halves based upon the rank order expression of this gene set (p = 0.06) (Fig. 3B). When tested in multivariate Cox proportional hazard models along with intrinsic subtype classification, MVD scores significantly contributed metastasis prediction information (p = 0.04).Fig. 3

Gene expression signatures add prognostic information to immunohistochemistry defined microvessel density scores. a Box-and-whisker plots are shown for average microvessel density scores for 70 human breast tumors. b–e Kaplan–Meier plots for relapse free survival and log-rank test p-values. For testing more than one variable c–e, tumors were independently ranked from low to high signature score and then the two groups were combined, which yielded groups of not necessarily equal number that were reflective of the biology of the tumor. The p-value in d and e test the tumors with low signature scores for both variables against all other tumors

We next aimed to elucidate why claudin-low tumors had the highest expression of the gene signatures that were designed to measure total endothelial quantity: Wallgard et al., vascular content, and activated endothelium non-proliferation component (Fig. 1; Supplemental Fig. 1). We reasoned that high vascular signature scores could either be attributed to the amount of vasculature present in these tumors, or due to the extent in which claudin-low breast cancer cells express endothelial genes. Since, the MVD scores suggest similar amounts of vasculature across the subtypes as assessed histologically (Fig. 3A), we hypothesized that claudin-low tumor cells themselves may express these vascular cell associated genes. Therefore, we identified vascular signature scores for the human breast cancer cell lines presented in Neve et al. [2, 27.] and found that the vascular content gene expression signature was most predominately expressed in claudin-low as compared to basal-like (p < 0.01) or luminal (p < 0.001) breast cancer cell lines (Fig. 4A). To further test the hypothesis that claudin-low cell lines have endothelial cell characteristics, we performed unsupervised hierarchical clustering using gene expression data from six breast cancer cell lines representative of different intrinsic subtypes and blood vessel endothelial cells (BECs). Interestingly, clustering with all available expressed genes (12,644) (Fig. 4B), or the vascular content gene signature (not shown), showed that claudin-low cell lines (i.e. MDA-MB-231 and SUM159) are transcriptomically more similar to BECs than they are to other breast epithelial cancer lines.Fig. 4

Gene expression and morphologic relatedness of endothelial cells and breast cancer cell lines. a Box-and-whisker plots of vascular signatures found in human breast cancer cell lines: BL; basal-like, CL; claudin-low, LUM; luminal. b Unsupervised hierarchical cluster dendogram of breast cancer cell lines and endothelial cells using all variably expressed genes (n = 12,644). c Picture of each cell line after 18 h of 3D culture (×5). d Pictures of the FAC sorted SUM149 cell line fractions after 18 h of 3D culture (×5). e Pictures of cocultures after 18 h of 3D culture (×5)

Given the genomic and morphologic similarities of claudin-low cell lines and endothelial cell lines, we next aimed to determine if tubular structures were formed by claudin-low cancer cells (Fig. 5A–F), as compared to luminal cancer cells (Fig. 5G–I), growing in vivo. To identify if any tube-like structures were able to functionally perfuse blood, mice were injected intravenously with Texas Red labeled dextran 5 min before euthanasia [15]. When subjected to pan-endothelial antibodies platelet/endothelial cell adhesion molecule (PECAM), von Willebrand factor (vWF), and the lymphatic vessel endothelial hyaluronan receptor 1 (LYVE1) simultaneously, the claudin-low tumors were found to have  extensive perfusion of dextran through paracellular spaces (Fig. 5B, E). This heightened vascular permeability was not observed in the luminal MCF-7 model (Fig. 5H). Serial frozen sections that utilized the cancer cell markers vimentin (Fig. 5C, F) or CK19 (Fig. 5I) confirm that the dextran freely diffused throughout and around the claudin-low tumors but was largely restricted to the vasculature in luminal tumors.Fig. 5

Identification of paracellular perfusion in claudin-low tumors. Texas Red Dextran (red) was injected into the circulation of mice bearing MDA231 (a–c), SUM159 (d–f), or MCF7 (g–i) cells grown as tumors in vivo. Serial sections for each tumor are shown. Hematoxylin and eosin staining (a, d, g); pan-endothelial antibodies vWF/PECAM/LYVE (green), DAPI (blue), and dextran (red) (b, e, h); vimentin or cytokeratin 19 cancer cell markers (green), DAPI (blue), and dextran (red) (c, f, i). Brackets denote extensive dextran perfusion in the absence of vasculature. All images are ×20