Infantile hemangiomas are the most commonly occurring tumors in childhood (affecting approximately 5-10% of children) and are composed of non-malignant, largely cutaneous vascular cells displaying a strongly over-proliferative phenotype. While infantile hemangiomas can occur throughout the body, these tumors exhibit a predilection for regions of facial embryological primordial fusion [1]. The growth of these benign tumors is divided into distinct stages which include proliferating, involuting, and involuted [2]. The proliferating stage occurs during the first year of life where the tumors often achieve their maximum size, generally ranging from the size of a pin-prick to 20 cm (average of 2 to 5 cm), with rare cases affecting a large portion of the infant’s body. When the child reaches 12-15 months of age, the tumors will cease growth and a spontaneous slow involution of the hemangioma proceeds. Tumors of most children are fully involuted by 5 to 10 years of age. Following involution, normal skin texture is restored in only 50% of children, leading to at least some degree of disfigurement in the former tumor area in at least half of the affected children. While 80-90% of hemangiomas are small and not dangerous, the remainders pose considerable risks to these children. For example, infantile hemangiomas that are alarming due to size, site, debilitating location (eyelid and orbit, airway, etc.), or that are life threatening (massive tumor, ulceration, hemorrhagic, or located on the heart or other vital organs) will require active therapeutic management. Complications occur in many children with these tumors causing functional or cosmetic issues, with the most common being infection followed by ulceration (21%) and bleeding (7%), but sometimes resulting in temporary or permanent blindness, airway obstruction, feeding difficulties, and hearing problems [3]. In very rare cases large infantile hemangiomas of the face may induce overgrowth of the facial skeleton or cartilage, resulting in gross disfigurement. Due to health risks and the psychological ramifications this disease has on pediatric patients, families, caregivers, and the extended community, the impact of this disease is disproportionate to the numbers affected.

Unfortunately, very little is known regarding the pathogenesis of this disease largely due to the absence of an analogous animal model that accurately recapitulates the disease progression and the exceptional difficulty in obtaining tumor samples from pediatric patients for a condition which most clinicians consider to be a self limiting disease. Familial clustering of infantile hemangiomas has been reported in a retrospective case-control study [4], suggesting a potential genetic contribution to this disease. Interestingly, females are 3 times more likely to have an infantile hemangioma than males, and females are 9 times more likely to display complications from these tumors [5]. Mutations in key angiogenic genes such as Vegfr2, Dusp5, and Tem8 have been reported in some infantile hemangiomas and these mutations may contribute to increased endothelial cell proliferation via constitutive stimulation of signaling pathways downstream of vascular endothelial growth factor (VEGF) [6-8].

Our lab has recently reported that endothelial cells purified from infantile hemangioma tumors are unique in their expression pattern relative to normal resident dermal endothelial cells [9]. The published analysis of this dataset focused primarily on aberrant angiogenic and tumor regulators that were differentially expressed between the two cell populations. However a statistically powered bioinformatics analysis of this dataset is lacking. In the current study, we utilized powerful bioinformatics software packages to statistically rank gene process networks differentially expressed between endothelial cells from the normal human dermis and infantile hemangioma tumors. We then performed network analysis to reveal known protein-protein interactions between a core set of 29 genes that were identified as significantly overrepresented these gene process networks. Finally, we used bioinformatics to identify putative transcriptional factors that may play a major role in regulating the expression of this core set of 29 genes.

Microarray dataset

Our meta-analysis was based on a study carried out by Stiles et al. [9] and deposited in Gene Expression Omnibus (GSE43742). In this study, the authors compared the gene expression profiles of Primary Cultures of Infantile Hemangiomas Endothelial Cells (HIHECs) to Human Dermal Micro Vascular Endothelial Cells (HDMVECs). Differential expression analysis showed that 125 genes were up regulated and 104 genes were down regulated in the HIHECs relative to the HDMVECs.

Gene ontology and pathway analysis

To understand the functions of the gene list, we utilized MetaCore software (Thompson Reuters) which provides a collection of high-level functions and utilities that each gene is putatively associated with. To provide functional annotation to each gene, we used MetaCore to identify statistically over-represented cellular processes in the dataset.

Protein-protein interaction network construction

To demonstrate the potential relationships among the genes in the most over-represented cellular processes, we matched the list of differentially expressed genes to the STRING database of known and predicted protein interactions derived from genomic context, highthroughput experiments, coexpression, and literature databases (http://string-db.org), and generated network maps for the genes involved in over-represented cellular processes.

Stiles et al. [9] previously used microarray technology to identify 229 genes that were differentially regulated between HIHECs and HDMVECs. These genes were involved in a number of cellular processes such as cell adhesion, cell cycle, and arachidonic acid production. This analysis identified and focused primarily on a number of angiogenic and tumor regulators that were aberrantly expressed in HIHECs. In the current study, we sought to utilize computational network analysis to statistically identify the key gene networks that are differentially expressed between HIHECs and HDMVECs and extrapolate putative biological processes that are deregulated in HIHECs.

We utilized MetaCore software for functional meta-analysis of the process networks most significantly represented by the 229 significantly expressed genes in the dataset. This analysis identified gene networks involved in cell adhesion and extracellular matrix remodeling as the most statistically represented pathways differentially regulated between the HIHECs and HDMVECs (Table 1). Generation of network maps for each of these processes identified a strong central theme indicative of collagen regulation including increases in the expression of COL5A1, COL5A2, COL4A6, and MMP2 and decreased expression of COL1A1, COL1A2, MMP1/13, and TIMP1, suggesting that HIHECs exhibit altered adhesive properties via changes in the expression of collagens and their direct regulators. A representative diagram of the most significantly represented process network (cell adhesion-cell matrix interactions) is shown in Figure 1.

Table 1: Seven most statistically significant process network maps generated from the 229 significantly expressed genes in the Stiles dataset.

Figure 1: Metacore generated process map of the most significantly represented process network (cell adhesion-cell matrix interactions) identified in our dataset. The ten genes in the far left column are upregulated in HIHECs relative to HDMVECs, while the five genes in the far right column are downregulated, respectively.

Combining all of the differentially expressed genes from the top seven most represented gene networks resulted in a reduction of the original 229 significantly expressed genes to a core 29 genes that appear to be central to the biological distinction between HIHECs and HDMVECs (Table 2 and Figure 2). We used the protein-protein interaction data provided by the STRING database to infer the functional interactions between each of the 29 core genes and generated an interaction network map depicting the interactions of these genes and their protein products (Figure 3). This data revealed that 71% of these core genes showed a positive interaction centered primarily on collagen function and proteolysis, further substantiating that cellular adhesions to the extracellular matrix are markedly deregulated in HIHECs.

Table 2: Table of the 29 core genes identified as overrepresented in our network analysis.

Figure 2: Hiearchical clustered heatmap of the normalized intensity values of the 29 core genes identified as overrepresented in our network analysis. (Red=high intensity value, green=low intensity value).

Figure 3: A protein-protein interaction map was generated using the STRING database to illustrate the known interactions between the 29 core genes identified as overrepresented in our network analysis.

We then used MetaCore to build a list of transcription factors that are statistically overrepresented as reported regulators of these core 29 genes. The top seven putative transcription factors that likely control this gene process network include SP1, RelA, CREB, ESR1, STAT3, ETS1, c-Jun, and C/EBP-beta (Table 3). We show network maps depicting how ESR1, STAT3, and C/EBP-beta (which have been implicated in the pathogenesis of infantile hemangiomas) have been reported to influence the expression of our core set of genes (Figure 4).

Table 3: Putative transcription factors which may control the core set of 29 genes based on statistical analysis on known regulation events.

Figure 4: Network map for the transcription factors ESR1, STAT3, and C/EBP-beta indicating known transcriptional regulations for the 29 core genes identified as overrepresented in our network analysis.

Our meta-analysis of data published by Stiles et al. [9], attempts to statistically identify the key gene networks that differentiate HIHECs from HDMVECs. We have used bioinformatics analysis to elucidate a core set of 29 genes whose involvement in biological processes such as cell adhesion and extracellular matrix degradation may help explain the aberrant signaling that drives the pathogenesis of infantile hemangiomas. Nearly three quarters of these core genes demonstrated known interactions, forming a tight network that largely centered on the regulation of collagen and extracellular matrix proteolysis.

The upregulation of COL5A1, COL5A2, COL4A6, and MMP2 and decreased expression of COL1A1, COL1A2, MMP1/13, and TIMP1 were the most interconnected hubs in our network analysis. The collagens COL5A1, COL5A2, COL4A6, COL1A1, and COL1A2 form the fibrillar acellular components of basement membranes and connective tissue [10]. These proteins play essential roles in both development and homeostasis across all tissues, and their aberrant expression is a hallmark of diseases such as cancer [11,12]. In addition to collagen, HIHECs displayed altered expression of FBN2, a fibrillin protein that associated with elastins in the cell matrix [13]. Matrix metalloproteinases such as MMP1/13 and MMP2 are responsible for proteolysis of the collagens, while the tissue inhibitor of metalloproteinases TIMP1 is a natural inhibitor of the MMPs [14] and altered expression of these proteins in HIHECs suggest strong deregulation of matrix cleavage. Other identified core genes involved in extracellular matrix adhesion include the hyaluronic acid receptors CD44 and LYVE1 [15,16], and CAV1, which has been shown to link integrin attachment sites to the Ras/ERK pathway [17]. The interplay and tight regulation of these proteins serves a central role in cell proliferation, migration, differentiation, apoptosis, angiogenesis, and host defenses [18]. Our data substantiates previous studies that have demonstrated altered extracellular remodeling in HIHECs via deregulated expression of collagens and MMPs [19-22].

Additionally, several genes involved in or affected by matrix proteolysis were included in the core gene set including F2RL1, ADAM19, UCHL1, SERPINE2, and SERPINB2. F2RL1 encodes a protease-activated thrombin receptor that is activated by trypsinmediated cleavage and has been shown to play a major role in vascular tone and hypotension via coupling to the nitric oxide pathway, as well as vessel inflammation and wound healing by stimulating the activation of immune cells and regulating capillary permeability [23-25]. ADAM19 encodes a membrane anchored metalloproteinase that participates in the proteolytic processing of beta-type neuregulin isoforms and plays a role in a variety of biological processes including cell-cell and cell-matrix interactions [26]. UCHL1 encodes a thiol protease belonging to the peptidase C12 family which hydrolyzes peptide bonds at the C-terminal glycine of ubiquitin and interacts with a number of NCAM cadherin proteins [27]. Both SERPINE2 and SERPINB2 encode inhibitors of serine proteases including thrombin, urokinase, plasmin, and trypsin [28], further emphasizing the role of extracellular proteolysis as a defining factor distinguishing HIHECs.

Our analysis identified putative transcription factors that may play a role in regulating this core gene network in HIHECs, including SP1, RelA, CREB, ESR1, STAT3, ETS1, c-Jun, and C/EBP-beta. While no reports in the literature have indicated that SP1, RelA, CREB, ETS1, or c-Jun are deregulated in infantile hemangiomas, C/EBP-beta, STAT3, and ESR1 have been shown to impinge on the pathology of these benign tumors. C/EBP-beta is one of the most important transcriptional regulators for collagen genes [29] and plays a major role in adipogenesis (a key process in the involution of infantile hemangioma tumors). Propranolol, the gold standard treatment for infantile hemangiomas, has been shown to enhance adipogenesis via increasing the expression levels of the pro-adipogenic transcription factor C/EBP-beta and other C/EBP family members [30]. STAT3 is an important marker for human embryonic stem cells (hESC) [31]. STAT3 along with Oct-4 has recently been shown to be expressed strongly in the endothelium of proliferating infantile hemangiomas [32], suggesting the presence of a primitive cellular origin for these tumors downstream of hESCs. ESR1, which encodes the estrogen receptor alpha protein, plays a predominant role in the vasculature by controlling cell survival, inflammation, oxidative stress, vasodilation, and neovascularization [33]. Patients with infantile hemangiomas demonstrate significantly higher estradiol (E2) levels compared to controls [34] and the combination of hypoxia and estrogen has been shown to synergistically enhance hemangioma endothelial cell proliferation [35]. We have previously published that the SWI/SNF chromatin remodeling subunit BAF57/SMARCE1 is significantly upregulated in HIHECs relative to HDMVECs [9], and this protein has been shown to be a critical regulator of estrogen receptor function in a number of cell types [36]. The estrogen link to infantile hemangiomas may help to explain why female infants are 3 times more likely to have an infantile hemangioma and 9 times more likely to show complications from these tumors.

In summary, our meta-analysis strongly indicates that changes in extracellular matrix proteolysis and cell adhesion are some of the most significantly defining factors influencing the pathology of infantile hemangiomas. Therapeutic interventions which aim to target these processes may show great efficacy in the treatment of patients afflicted with this tumor.

Funding for this work was kindly provided by TTUHSC startup funds to BAB and SARP summer medical student research scholarships to DT and KV.