Microarray and ChIP-seq data analysis revealed changes in p53-mediated transcriptional regulation in Nutlin-3-treated U2OS cells

Integrative analysis of chromatin immunoprecipitation-sequencing (ChIP-seq) data and microarray data was performed to illustrate the effect of Nutlin-3 on promoter selectivity and transcriptional regulation by the tumor suppressor p53 in U2OS human osteosarcoma cells. Raw data (accession number, GSE46642) were downloaded from Gene Expression Omnibus. Differential analyses were performed using package limma of R software. Gene ontology enrichment and Kyoto Encyclopedia of Genes and Genomes pathway enrichment analyses were performed for the differentially expressed genes (DEGs) using the Database for Annotation, Visualization and Integration Discovery. Integrative analysis of ChIP-seq data and microarray data were confirmed with ChIP-Array. A total of 565 DEGs were identified, including 373 upregulated genes and 192 downregulated genes. Genes involved in the p53 signaling pathway, cell cycle, DNA replication, cytokine-cytokine receptor interaction and melanoma were markedly over-represented in the DEGs. A total of 39 DEGs were directly regulated by p53 and two were the transcription factors (TFs), E2F2 and HOXA1. E2F2 regulated 25 DEGs, while HOXA1 regulated one DEG. The cell cycle, p53 signaling pathway, melanoma and pathways involved in cancer were enriched in the direct and indirect target genes. Changes in the p53-binding pattern induced by Nutlin-3 were described in the present study, which may advance the understanding of the regulatory network of p53 in osteosarcoma and aid in the development of novel therapies.


Introduction
Gene transcription is regulated by dynamic interactions between cis-regulatory elements and regulatory proteins, including transcription factors (TFs). Tumor protein p53 is an important TF involved in various cellular processes, including growth arrest, senescence and apoptosis (1)(2)(3). Following cellular stress, stabilized p53 translocates into the nucleus and subsequently binds to the consensus sequence motif to regulate the expression of hundreds of genes.
p53 is critical in tumor suppression and loss of p53 function is required for cancer progression. Mutational inactivation of p53 is detected in >50% of human cancer types (4). A number of downstream proteins of p53 have been identified (5)(6)(7). Nevertheless, several of the factors expected to affect p53-induced changes in gene expression are poorly understood, including the impact of different stresses that can induce p53. Genome-wide studies may provide an improved understanding of its transcriptional regulatory functions in certain types of cancer (8)(9)(10), including osteosarcoma.
Osteosarcoma is the eighth most common type of childhood cancer and is also the most common histological form of primary bone cancer (11). The mortality rates for osteosarcoma have been declining by ~1.3% annually (12). The overall 5-year survival rate for osteosarcoma is ~68% (12). Future studies are required to fully disclose the molecular mechanisms and advance therapeutic development.
In the present study, human U2OS osteosarcoma cells, expressing wild-type p53, were used to investigate the effect of treatment with Nutlin-3 (a non-genotoxic activator of p53) on p53 binding genes. Different from a previous study by Menendez et al (13), a stricter threshold [|log 2 fold change (FC)|>1 and false discovery rate (FDR) <0.05 vs. FC>2 and FDR ≤0.1] was used to select the differentially expressed genes (DEGs) and to construct the regulatory association between p53 and its target genes.

Materials and methods
Raw data. The raw data (accession number, GSE46642) were downloaded from Gene Expression Omnibus (http://www.
Pre-treatment and differential analysis. The microarray data were read using the package, affy (14), on the software R (http://www.r-project.org/). Following background correction and normalization with a Robust Multi-array Analysis (RMA) method in R affy, the gene expression levels were determined. Differential analysis was performed using the package, linear models for microarray data (limma) (15), on the software R.
Multiple-testing correction was performed using the Bayes method (implemented in the 'limma' R package). The following threshold was set for the screening of the DEGs: |log 2 FC|>1 and FDR<0.05.

Integrative analysis of microarray data and ChIP-seq data.
ChIP-Array (http://jjwanglab.org/chip-array) is an online tool developed for integrative analysis of microarray data and ChIP-seq data (16). It identifies the indirect target, Z, by identifying an intermediate transcription factor (TF), Y, which is a putative regulator of Z and a target of X. The putative regulator of Z is identified by scanning all promoters in the genome with position weight matrix (PWMs) of all Ys from three publicly accessible databases [JASPAR (http://jaspar.genereg. net), UniPROBE (http://uniprobe.org) and TRANSFAC (http://www.gene-regulation.com/pub/databases.html) derived transcription factor binding site database from University of California, Santa Cruz genome browser] (16).
In the present study, the parameters were set as follows: Promoter range, -500~+100; TF database, UniPROBE; PWM scan P-value, 10 -5 ; and conservation filtering P-value, 0.001. Finally, a gene regulatory network was obtained for p53, including its direct and indirect target genes.

Functional enrichment analysis. Gene Ontology (GO) and
Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses were performed for the DEGs using the Database for Annotation, Visualization and Integration Discovery (http://david.abcc.ncifcrf.gov/) online tools (17). P<0.05 was considered to indicate a statistically significant difference and was set as the cut-off.

Results
Differentially expressed genes. Gene expression data prior to and following normalization with the RMA method are demonstrated in Fig. 1. A good performance of normalization was achieved.
A total of 565 DEGs were identified, including 373 upregulated genes and 192 downregulated genes. Clustering and a heat-map of the expression values for DEGs are shown in Fig. 2. The Nutlin-3 treated U2OS samples were well distinguished from the control samples, suggesting the reliability of the DEGs. Table I. Kyoto Encyclopedia of Genes and Genomes pathways enriched in the differentially expressed genes. Pathway  Table I. The p53 signaling pathway, cell cycle, DNA replication, cytokine-cytokine receptor interaction and melanoma were significantly over-represented in the DEGs.
Transcriptional regulatory network of p53. Integrative analysis of ChIP-seq data and microarray data was performed using the ChIP-Array online tool. A total of 39 DEGs were directly regulated by p53, and two of them were TFs: E2F transcription factor 2 (E2F2) and homeobox A1 (HOXA1). E2F2 regulated 25 DEGs and HOXA1 regulated one DEG (Fig. 4).
Functional enrichment analysis result of the target genes. GO enrichment analysis was performed for the direct and indirect target genes of p53 (Fig. 5). Cell cycle and cell-cell signaling were included in the list.
The KEGG pathway enriched in all the target genes of p53 were also disclosed (Table II), including cell cycle, p53 signaling pathway, melanoma and pathways in cancer.

Discussion
In the present study, a total of 565 DEGs were identified in Nutlin-3-treated U2OS cells compared with the control samples. Of these DEGs, 373 were upregulated genes and 192 were downregulated genes. Functional enrichment analysis revealed that the p53 signaling pathway, cell cycle and DNA replication were significantly over-represented in the DEGs. This result suggested the importance of p53 in osteosarcoma. p53 functions as a cell cycle control protein in osteosarcoma (18) and the presence of p53 mutations in human osteosarcoma is correlated with high levels of genomic instability (19), confirming the critical importance of p53 in response to stresses, including DNA damage. Berman et al (20) reported that metastatic osteosarcoma is induced by the inactivation of Rb and p53 (20). The comparative analysis of gene expression profiles between Nutlin-3-treated U2OS cells and controls further described the critical importance of p53 in osteosarcoma. Notably, p53 gene therapy of human osteosarcoma is also suggested and has been previously investigated (21).
In conclusion, differential expression of several direct and indirect target genes of p53 was observed following treatment with Nutlin-3. These findings not only advanced the understanding regarding the importance of p53 in osteosarcoma, but also provided clues for future development of therapeutic strategies.