Reusable, extensible, and modifiable R scripts and Kepler workflows for comprehensive single set ChIP-seq analysis.

BACKGROUND: There has been an enormous expansion of use of chromatin immunoprecipitation followed by sequencing (ChIP-seq) technologies. Analysis of large-scale ChIP-seq datasets involves a complex series of steps and production of several specialized graphical outputs. A number of systems have emphasized custom development of ChIP-seq pipelines. These systems are primarily based on custom programming of a single, complex pipeline or supply libraries of modules and do not produce the full range of outputs commonly produced for ChIP-seq datasets. It is desirable to have more comprehensive pipelines, in particular ones addressing common metadata tasks, such as pathway analysis, and pipelines producing standard complex graphical outputs. It is advantageous if these are highly modular systems, available as both turnkey pipelines and individual modules, that are easily comprehensible, modifiable and extensible to allow rapid alteration in response to new analysis developments in this growing area. Furthermore, it is advantageous if these pipelines allow data provenance tracking. RESULTS: We present a set of 20 ChIP-seq analysis software modules implemented in the Kepler workflow system; most (18/20) were also implemented as standalone, fully functional R scripts. The set consists of four full turnkey pipelines and 16 component modules. The turnkey pipelines in Kepler allow data provenance tracking. Implementation emphasized use of common R packages and widely-used external tools (e.g., MACS for peak finding), along with custom programming. This software presents comprehensive solutions and easily repurposed code blocks for ChIP-seq analysis and pipeline creation. Tasks include mapping raw reads, peakfinding via MACS, summary statistics, peak location statistics, summary plots centered on the transcription start site (TSS), gene ontology, pathway analysis, and de novo motif finding, among others. CONCLUSIONS: These pipelines range from those performing a single task to those performing full analyses of ChIP-seq data. The pipelines are supplied as both Kepler workflows, which allow data provenance tracking, and, in the majority of cases, as standalone R scripts. These pipelines are designed for ease of modification and repurposing.

An inhaled dose of budesonide induces genes involved in transcription and signaling in the human airways: enhancement of anti- and proinflammatory effector genes.

Although inhaled glucocorticoids, or corticosteroids (ICS), are generally effective in asthma, understanding their anti-inflammatory actions in vivo remains incomplete. To characterize glucocorticoid-induced modulation of gene expression in the human airways, we performed a randomized placebo-controlled crossover study in healthy male volunteers. Six hours after placebo or budesonide inhalation, whole blood, bronchial brushings, and endobronchial biopsies were collected. Microarray analysis of biopsy RNA, using stringent (≥2-fold, 5% false discovery rate) or less stringent (≥1.25-fold, P ≤ 0.05) criteria, identified 46 and 588 budesonide-induced genes, respectively. Approximately two third of these genes are transcriptional regulators (KLF9, PER1, TSC22D3, ZBTB16), receptors (CD163, CNR1, CXCR4, LIFR, TLR2), or signaling genes (DUSP1, NFKBIA, RGS1, RGS2, ZFP36). Listed genes were qPCR verified. Expression of anti-inflammatory and other potentially beneficial genes is therefore confirmed and consistent with gene ontology (GO) terms for negative regulation of transcription and gene expression. However, GO terms for transcription, signaling, metabolism, proliferation, inflammatory responses, and cell movement were also associated with the budesonide-induced genes. The most enriched functional cluster indicates positive regulation of proliferation, locomotion, movement, and migration. Moreover, comparison with the budesonide-induced expression profile in primary human airway epithelial cells shows considerable cell type specificity. In conclusion, increased expression of multiple genes, including the transcriptional repressor, ZBTB16, that reduce inflammatory signaling and gene expression, occurs in the airways and blood and may contribute to the therapeutic efficacy of ICS. This provides a previously lacking insight into the in vivo effects of ICS and should promote strategies to improve glucocorticoid efficacy in inflammatory diseases.

Differences among brain tumor stem cell types and fetal neural stem cells in focal regions of histone modifications and DNA methylation, broad regions of modifications, and bivalent promoters.

BACKGROUND: Aberrational epigenetic marks are believed to play a major role in establishing the abnormal features of cancer cells. Rational use and development of drugs aimed at epigenetic processes requires an understanding of the range, extent, and roles of epigenetic reprogramming in cancer cells. Using ChIP-chip and MeDIP-chip approaches, we localized well-established and prevalent epigenetic marks (H3K27me3, H3K4me3, H3K9me3, DNA methylation) on a genome scale in several lines of putative glioma stem cells (brain tumor stem cells, BTSCs) and, for comparison, normal human fetal neural stem cells (fNSCs). RESULTS: We determined a substantial "core" set of promoters possessing each mark in every surveyed BTSC cell type, which largely overlapped the corresponding fNSC sets. However, there was substantial diversity among cell types in mark localization. We observed large differences among cell types in total number of H3K9me3+ positive promoters and peaks and in broad modifications (defined as >50 kb peak length) for H3K27me3 and, to a lesser extent, H3K9me3. We verified that a change in a broad modification affected gene expression of CACNG7. We detected large numbers of bivalent promoters, but most bivalent promoters did not display direct overlap of contrasting epigenetic marks, but rather occupied nearby regions of the proximal promoter. There were significant differences in the sets of promoters bearing bivalent marks in the different cell types and few consistent differences between fNSCs and BTSCs. CONCLUSIONS: Overall, our "core set" data establishes sets of potential therapeutic targets, but the diversity in sets of sites and broad modifications among cell types underscores the need to carefully consider BTSC subtype variation in epigenetic therapy. Our results point toward substantial differences among cell types in the activity of the production/maintenance systems for H3K9me3 and for broad regions of modification (H3K27me3 or H3K9me3). Finally, the unexpected diversity in bivalent promoter sets among these multipotent cells indicates that bivalent promoters may play complex roles in the overall biology of these cells. These results provide key information for forming the basis for future rational drug therapy aimed at epigenetic processes in these cells.

Mutation in folate metabolism causes epigenetic instability and transgenerational effects on development.

The importance of maternal folate consumption for normal development is well established, yet the molecular mechanism linking folate metabolism to development remains poorly understood. The enzyme methionine synthase reductase (Mtrr) is necessary for utilization of methyl groups from the folate cycle. We found that a hypomorphic mutation of the mouse Mtrr gene results in intrauterine growth restriction, developmental delay, and congenital malformations, including neural tube, heart, and placental defects. Importantly, these defects were dependent upon the Mtrr genotypes of the maternal grandparents. Furthermore, we observed widespread epigenetic instability associated with altered gene expression in the placentas of wild-type grandprogeny of Mtrr-deficient maternal grandparents. Embryo transfer experiments revealed that Mtrr deficiency in mice lead to two distinct, separable phenotypes: adverse effects on their wild-type daughters' uterine environment, leading to growth defects in wild-type grandprogeny, and the appearance of congenital malformations independent of maternal environment that persist for five generations, likely through transgenerational epigenetic inheritance.

Workflows for microarray data processing in the Kepler environment.

BACKGROUND: Microarray data analysis has been the subject of extensive and ongoing pipeline development due to its complexity, the availability of several options at each analysis step, and the development of new analysis demands, including integration with new data sources. Bioinformatics pipelines are usually custom built for different applications, making them typically difficult to modify, extend and repurpose. Scientific workflow systems are intended to address these issues by providing general-purpose frameworks in which to develop and execute such pipelines. The Kepler workflow environment is a well-established system under continual development that is employed in several areas of scientific research. Kepler provides a flexible graphical interface, featuring clear display of parameter values, for design and modification of workflows. It has capabilities for developing novel computational components in the R, Python, and Java programming languages, all of which are widely used for bioinformatics algorithm development, along with capabilities for invoking external applications and using web services. RESULTS: We developed a series of fully functional bioinformatics pipelines addressing common tasks in microarray processing in the Kepler workflow environment. These pipelines consist of a set of tools for GFF file processing of NimbleGen chromatin immunoprecipitation on microarray (ChIP-chip) datasets and more comprehensive workflows for Affymetrix gene expression microarray bioinformatics and basic primer design for PCR experiments, which are often used to validate microarray results. Although functional in themselves, these workflows can be easily customized, extended, or repurposed to match the needs of specific projects and are designed to be a toolkit and starting point for specific applications. These workflows illustrate a workflow programming paradigm focusing on local resources (programs and data) and therefore are close to traditional shell scripting or R/BioConductor scripting approaches to pipeline design. Finally, we suggest that microarray data processing task workflows may provide a basis for future example-based comparison of different workflow systems. CONCLUSIONS: We provide a set of tools and complete workflows for microarray data analysis in the Kepler environment, which has the advantages of offering graphical, clear display of conceptual steps and parameters and the ability to easily integrate other resources such as remote data and web services.

Anesthetics discriminate between tonic and phasic gamma-aminobutyric acid receptors on hippocampal CA1 neurons.

BACKGROUND: Anesthesia is produced by a depression of neuronal signaling in the central nervous system (CNS); however, the mechanism(s) of action underlying this depression remain unclear. Recent studies have indicated that anesthetics can enhance inhibition of CNS neurons by increasing current flow through tonic gamma-aminobutyric acid (GABA(A)) receptor gated chloride channels in their membranes. Enhanced tonic inhibition would contribute to CNS depression produced by anesthetics, but it remains to be determined to what extent anesthetic actions at these receptors contribute to CNS depression. In the present study, we compared and contrasted the involvement of tonic versus synaptic GABA(A) receptors in the functional depression of CNS neurons produced by isoflurane and thiopental. METHODS: In rat hippocampal slices, whole cell patch clamp recordings were used to study anesthetic effects on CA1 neuron intrinsic excitability, and population spike recordings were used to investigate effects on synaptically evoked discharge. These responses were chosen to test whether anesthetic effects on GABA receptors alter single neuron discharge and/or circuit level synaptic functioning. Phasic (synaptic) GABA receptors were selectively blocked using the GABA(A) antagonist gabazine and tonic responses were blocked using the chloride channel blocker picrotoxin. RESULTS: Clinically relevant and equi-effective concentrations of thiopental and isoflurane depressed CA1 neuron synaptically evoked discharge. This depression was partially reversed by blocking synaptic GABA(A) receptors with gabazine (20 microM). The thiopental-induced depression was reversed by approximately 60%, but the isoflurane-induced depression was reversed by only approximately 20%. Blocking tonic GABA(A) receptors with the addition of 100 microM picrotoxin produced an additional 40% reversal of the thiopental-induced depression, but no additional reversal was seen for isoflurane-depressed responses. In response to direct DC current injection, CA1 neuron discharge was depressed by thiopental and membrane conductance was increased. Both of these effects were reversed by picrotoxin, but not by gabazine. Isoflurane, in contrast, neither depressed current-evoked discharge, nor altered the membrane conductance of CA1 neurons. CONCLUSIONS: These results indicate that general anesthetics discriminate between synaptic and tonic GABA(A) receptors. Effects on both phasic and tonic receptors combined to depress circuit responses produced by thiopental, whereas only effects on synaptic GABA receptors appeared to play an important role for isoflurane. Together with the other known sites of action for these two anesthetics, our results support a multisite, agent-specific mechanism for anesthetic actions.

Analysis of the mechanisms mediating tumor-specific changes in gene expression in human liver tumors.

There is widespread interest in efficient characterization of differences between tumor and normal samples. Here, we show an effective methodology for genome-scale characterization of tumors. Using matched normal and tumor samples from liver cancer patients, as well as non-cancer-related normal liver tissue, we first determined changes in gene expression as monitored on RNA expression arrays. We identified several hundred mRNAs that were consistently changed in the tumor samples. To characterize the mechanisms responsible for creation of the tumor-specific transcriptome, we performed chromatin immunoprecipitation on microarray experiments to assay binding of RNA polymerase II, H3me3K27, and H3me3K9 and DNA methylation in 25,000 promoter regions. These experiments identified changes in active and silenced regions of the genome in the tumor cells. Finally, we used a "virtual comparative genomic hybridization" method to identify copy number alterations in the tumor samples. Through comparison of RNA polymerase II binding, chromatin structure, DNA methylation, and copy number changes, we suggest that the major contributor to creation of the liver tumor transcriptome was changes in gene copy number.

Systematic evaluation of variability in ChIP-chip experiments using predefined DNA targets.

The most widely used method for detecting genome-wide protein-DNA interactions is chromatin immunoprecipitation on tiling microarrays, commonly known as ChIP-chip. Here, we conducted the first objective analysis of tiling array platforms, amplification procedures, and signal detection algorithms in a simulated ChIP-chip experiment. Mixtures of human genomic DNA and "spike-ins" comprised of nearly 100 human sequences at various concentrations were hybridized to four tiling array platforms by eight independent groups. Blind to the number of spike-ins, their locations, and the range of concentrations, each group made predictions of the spike-in locations. We found that microarray platform choice is not the primary determinant of overall performance. In fact, variation in performance between labs, protocols, and algorithms within the same array platform was greater than the variation in performance between array platforms. However, each array platform had unique performance characteristics that varied with tiling resolution and the number of replicates, which have implications for cost versus detection power. Long oligonucleotide arrays were slightly more sensitive at detecting very low enrichment. On all platforms, simple sequence repeats and genome redundancy tended to result in false positives. LM-PCR and WGA, the most popular sample amplification techniques, reproduced relative enrichment levels with high fidelity. Performance among signal detection algorithms was heavily dependent on array platform. The spike-in DNA samples and the data presented here provide a stable benchmark against which future ChIP platforms, protocol improvements, and analysis methods can be evaluated.

Integrated epigenomic analyses of neuronal MeCP2 reveal a role for long-range interaction with active genes.

Mutations in MECP2 cause the autism-spectrum disorder Rett syndrome. MeCP2 is predicted to bind to methylated promoters and silence transcription. However, the first large-scale mapping of neuronal MeCP2-binding sites on 26.3 Mb of imprinted and nonimprinted loci revealed that 59% of MeCP2-binding sites are outside of genes and that only 6% are in CpG islands. Integrated genome-wide promoter analysis of MeCP2 binding, CpG methylation, and gene expression revealed that 63% of MeCP2-bound promoters are actively expressed and that only 6% are highly methylated. These results indicate that the primary function of MeCP2 is not the silencing of methylated promoters.

A comprehensive ChIP-chip analysis of E2F1, E2F4, and E2F6 in normal and tumor cells reveals interchangeable roles of E2F family members.

Using ChIP-chip assays (employing ENCODE arrays and core promoter arrays), we examined the binding patterns of three members of the E2F family in five cell types. We determined that most E2F1, E2F4, and E2F6 binding sites are located within 2 kb of a transcription start site, in both normal and tumor cells. In fact, the majority of promoters that are active (as defined by TAF1 or POLR2A binding) in GM06990 B lymphocytes and Ntera2 carcinoma cells were also bound by an E2F. This very close relationship between E2F binding sites and binding sites for general transcription factors in both normal and tumor cells suggests that a chromatin-bound E2F may be a signpost for active transcription initiation complexes. In general, we found that several E2Fs bind to a given promoter and that there is only modest cell type specificity of the E2F family. Thus, it is difficult to assess the role of any particular E2F in transcriptional regulation, due to extreme redundancy of target promoters. However, Ntera2 carcinoma cells were exceptional in that a large set of promoters were bound by E2F6, but not by E2F1 or E2F4. It has been proposed that E2F6 contributes to gene silencing by recruiting enzymes involved in methylating histone H3. To test this hypothesis, we created Ntera2 cell lines harboring shRNAs to E2F6. We found that reduction of E2F6 only induced minimal alteration of the transcriptome of Ntera2 transcriptome. Our results support the concept of functional redundancy in the E2F family and suggest that E2F6 is not critical for histone methylation.

Identification of genes directly regulated by the oncogene ZNF217 using chromatin immunoprecipitation (ChIP)-chip assays.

It has been proposed that ZNF217, which is amplified at 20q13 in various tumors, plays a key role during neoplastic transformation. ZNF217 has been purified in complexes that contain repressor proteins such as CtBP2, suggesting that it acts as a transcriptional repressor. However, the function of ZNF217 has not been well characterized due to a lack of known target genes. Using a global chromatin immunoprecipitation (ChIP)-chip approach, we identified thousands of ZNF217 binding sites in three tumor cell lines (MCF7, SW480, and Ntera2). Further analysis of ZNF217 in Ntera2 cells showed that many promoters are bound by ZNF217 and CtBP2 and that a subset of these promoters are activated upon removal of ZNF217. Thus, our in vivo studies corroborate the in vitro biochemical analyses of ZNF217-containing complexes and support the hypothesis that ZNF217 functions as a transcriptional repressor. Gene ontology analysis showed that ZNF217 targets in Ntera2 cells are involved in organ development, suggesting that one function of ZNF217 may be to repress differentiation. Accordingly we show that differentiation of Ntera2 cells with retinoic acid led to down-regulation of ZNF217. Our identification of thousands of ZNF217 target genes will enable further studies of the consequences of aberrant expression of ZNF217 during neoplastic transformation.

Unbiased location analysis of E2F1-binding sites suggests a widespread role for E2F1 in the human genome.

The E2F family of transcription factors regulates basic cellular processes. Here, we take an unbiased approach towards identifying E2F1 target genes by examining localization of E2F1-binding sites using high-density oligonucleotide tiling arrays. To begin, we developed a statistically-based methodology for analysis of ChIP-chip data obtained from arrays that represent 30 Mb of the human genome. Using this methodology, we identified regions bound by E2F1, MYC, and RNA Polymerase II (POLR2A). We found a large number of binding sites for all three factors; extrapolation suggests there may be approximately 20,000-30,000 E2F1- and MYC-binding sites and approximately 12,000-17,000 active promoters in HeLa cells. In contrast to our results for MYC, we find that the majority of E2F1-binding sites (>80%) are located in core promoters and that 50% of the sites overlap transcription starts. Only a small fraction of E2F1 sites possess the canonical binding motif. Surprisingly, we found that approximately 30% of genes in the 30-Mb region possessed an E2F1 binding site in a core promoter and E2F1 was bound near to 83% of POLR2A-bound sites. To determine if these results were representative of the entire human genome, we performed ChIP-chip analyses of approximately 24,000 promoters and confirmed that greater than 20% of the promoters were bound by E2F1. Our results suggest that E2F1 is recruited to promoters via a method distinct from recognition of the known consensus site and point toward a new understanding of E2F1 as a factor that contributes to the regulation of a large fraction of human genes.

N-type and L-type calcium channels mediate glycinergic synaptic inputs to retinal ganglion cells of tiger salamanders.

Synaptically localized calcium channels shape the timecourse of synaptic release, are a prominent site for neuromodulation, and have been implicated in genetic disease. In retina, it is well established that L-type calcium channels play a major role in mediating release of glutamate from the photoreceptors and bipolar cells. However, little is known about which calcium channels are coupled to synaptic exocytosis of glycine, which is primarily released by amacrine cells. A recent report indicates that glycine release from spiking AII amacrine cells relies exclusively upon L-type calcium channels. To identify calcium channel types controlling neurotransmitter release from the population of glycinergic neurons that drive retinal ganglion cells, we recorded electrical and potassium evoked inhibitory synaptic currents (IPSCs) from these postsynaptic neurons in retinal slices from tiger salamanders. The L-channel antagonist nifedipine strongly inhibited release and FPL64176, an L-channel agonist, greatly enhanced it, indicating a significant role for L-channels. omega-Conotoxin MVIIC, an N/P/Q-channel antagonist, strongly inhibited release, indicating an important role for non-L channels. While the P/Q-channel blocker omega-Aga IVA produced only small effects, the N-channel blocker omega-conotoxin GVIA strongly inhibited release. Hence, N-type and L-type calcium channels appear to play major roles, overall, in mediating synaptic release of glycine onto retinal ganglion cells.

Major role for tonic GABAA conductances in anesthetic suppression of intrinsic neuronal excitability.

Anesthetics appear to produce neurodepression by altering synaptic transmission and/or intrinsic neuronal excitability. Propofol, a widely used anesthetic, has proposed effects on many targets, ranging from sodium channels to GABA(A) inhibition. We examined effects of propofol on the intrinsic excitability of hippocampal CA1 neurons (primarily interneurons) recorded from adult rat brain slices. Propofol strongly depressed action potential production induced by DC injection, synaptic stimulation, or high-potassium solutions. Propofol-induced depression of intrinsic excitability was completely reversed by bicuculline and picrotoxin but was strychnine-insensitive, implicating GABA(A) but not glycine receptors. Propofol strongly enhanced inhibitory postsynaptic currents (IPSCs) and induced a tonic GABA(A)-mediated current. We pharmacologically differentiated tonic and phasic (synaptic) GABA(A)-mediated inhibition using the GABA(A) receptor antagonist SR95531 (gabazine). Gabazine (20 microM) completely blocked both evoked and spontaneous IPSCs but failed to block the propofol-induced depression of intrinsic excitability, implicating tonic, but not phasic, GABA(A) inhibition. Glutamatergic synaptic responses were not altered by propofol (< or =30 microM). Similar results were found in both interneurons and pyramidal cells and with the chemically unrelated anesthetic thiopental. These results suggest that suppression of CA1 neuron intrinsic excitability, by these anesthetics, is largely due to activation of tonic GABA(A) conductances; although other sites of action may play important roles in affecting synaptic transmission, which also can produce strong neurodepression. We propose that for some anesthetics, suppression of intrinsic excitability, mediated by tonic GABA(A) conductances, operates in conjunction with effects on synaptic transmission, mediated by other mechanisms, to depress hippocampal function during anesthesia.