Efficiently identifying genome-wide changes with next-generation sequencing data

Comprehensive assessment of differential ChIP-seq tools guides optimal algorithm selection

Background

The analysis of chromatin binding patterns of proteins in different biological states is a main application of chromatin immunoprecipitation followed by sequencing (ChIP-seq). A large number of algorithms and computational tools for quantitative comparison of ChIP-seq datasets exist, but their performance is strongly dependent on the parameters of the biological system under investigation. Thus, a systematic assessment of available computational tools for differential ChIP-seq analysis is required to guide the optimal selection of analysis tools based on the present biological scenario.

Results

We created standardized reference datasets by in silico simulation and sub-sampling of genuine ChIP-seq data to represent different biological scenarios and binding profiles. Using these data, we evaluated the performance of 33 computational tools and approaches for differential ChIP-seq analysis. Tool performance was strongly dependent on peak size and shape as well as on the scenario of biological regulation.

Conclusions

Our analysis provides unbiased guidelines for the optimized choice of software tools in differential ChIP-seq analysis.

Supplementary Information

The online version contains supplementary material available at 10.1186/s13059-022-02686-y.

AnnoMiner is a new web-tool to integrate epigenetics, transcription factor occupancy and transcriptomics data to predict transcriptional regulators

Gene expression regulation requires precise transcriptional programs, led by transcription factors in combination with epigenetic events. Recent advances in epigenomic and transcriptomic techniques provided insight into different gene regulation mechanisms. However, to date it remains challenging to understand how combinations of transcription factors together with epigenetic events control cell-type specific gene expression. We have developed the AnnoMiner web-server, an innovative and flexible tool to annotate and integrate epigenetic, and transcription factor occupancy data. First, AnnoMiner annotates user-provided peaks with gene features. Second, AnnoMiner can integrate genome binding data from two different transcriptional regulators together with gene features. Third, AnnoMiner offers to explore the transcriptional deregulation of genes nearby, or within a specified genomic region surrounding a user-provided peak. AnnoMiner’s fourth function performs transcription factor or histone modification enrichment analysis for user-provided gene lists by utilizing hundreds of public, high-quality datasets from ENCODE for the model organisms human, mouse, Drosophila and C. elegans. Thus, AnnoMiner can predict transcriptional regulators for a studied process without the strict need for chromatin data from the same process. We compared AnnoMiner to existing tools and experimentally validated several transcriptional regulators predicted by AnnoMiner to indeed contribute to muscle morphogenesis in Drosophila. AnnoMiner is freely available at http://chimborazo.ibdm.univ-mrs.fr/AnnoMiner/.

Bioinformatic tools for DNA methylation and histone modification: A survey

Epigenetic inheritance occurs due to different mechanisms such as chromatin and histone modifications, DNA methylation and processes mediated by non-coding RNAs. It leads to changes in gene expressions and the emergence of new traits in different organisms in many diseases such as cancer. Recent advances in experimental methods led to the identification of epigenetic target sites in various organisms. Computational approaches have enabled us to analyze mass data produced by these methods. Next-generation sequencing (NGS) methods have been broadly used to identify these target sites and their patterns. By using these patterns, the emergence of diseases could be prognosticated. In this study, target site prediction tools for two major epigenetic mechanisms comprising histone modification and DNA methylation are reviewed. Publicly accessible databases are reviewed as well. Some suggestions regarding the state-of-the-art methods and databases have been made, including examining patterns of epigenetic changes that are important in epigenotypes detection.

Global transcriptomic profiling of microcystin-LR or -RR treated hepatocytes (HepaRG)

The canonical mode of action (MOA) of microcystins (MC) is the inhibition of protein phosphatases, but complete characterization of toxicity pathways is lacking. The existence of over 200 MC congeners complicates risk estimates worldwide. This work employed RNA-seq to provide an unbiased and comprehensive characterization of cellular targets and impacted cellular processes of hepatocytes exposed to either MC-LR or MC-RR congeners. The human hepatocyte cell line, HepaRG, was treated with three concentrations of MC-LR or -RR for 2 h. Significant reduction in cell survival was observed in LR1000 and LR100 treatments whereas no acute toxicity was observed in any MR-RR treatment. RNA-seq was performed on all treatments of MC-LR and -RR. Differentially expressed genes and pathways associated with oxidative and endoplasmic reticulum (ER) stress, and the unfolded protein response (UPR) were highly enriched by both congeners as were inflammatory pathways. Genes associated with both apoptotic and inflammatory pathways were enriched in LR1000. We present a model of MC toxicity that immediately causes oxidative stress and leads to ER stress and the activation of the UPR. Differential activation of the three arms of the UPR and the kinetics of JNK activation ultimately determine whether cell survival or apoptosis is favored. Extracellular exosomes were enrichment of by both congeners, suggesting a previously unidentified mechanism for MC-dependent extracellular signaling. The complement system was enriched only in MC-RR treatments, suggesting congener-specific differences in cellular effects. This study provided an unbiased snapshot of the early systemic hepatocyte response to MC-LR and MC-RR congeners and may explain differences in toxicity among MC congeners.

•

Microcystin-LR and microcystin-RR have similar transcriptional responses.

•

Genes associated with oxidative stress and the unfolded protein response were enriched by congeners.

•

Genes associated with extracellular exosomes were enriched, suggesting a potential new mechanism for cell signaling.

•

Complement associated genes were strongly enriched only by microcystin-RR.

•

Identified a potential molecular mechanism underlying the cellular fate of hepatocyte.

Coordinated regulation of Rel expression by MAP3K4, CBP, and HDAC6 controls phenotypic switching

Coordinated gene expression is required for phenotypic switching between epithelial and mesenchymal phenotypes during normal development and in disease states. Trophoblast stem (TS) cells undergo epithelial-mesenchymal transition (EMT) during implantation and placentation. Mechanisms coordinating gene expression during these processes are poorly understood. We have previously demonstrated that MAP3K4-regulated chromatin modifiers CBP and HDAC6 each regulate thousands of genes during EMT in TS cells. Here we show that CBP and HDAC6 coordinate expression of only 183 genes predicted to be critical regulators of phenotypic switching. The highest-ranking co-regulated gene is the NF-κB family member Rel. Although NF-κB is primarily regulated post-transcriptionally, CBP and HDAC6 control Rel transcript levels by binding Rel regulatory regions and controlling histone acetylation. REL re-expression in mesenchymal-like TS cells induces a mesenchymal-epithelial transition. Importantly, REL forms a feedback loop, blocking HDAC6 expression and nuclear localization. Together, our work defines a developmental program coordinating phenotypic switching.

Noha Shendy et al. study the role of CBP and HDAC6 in phenotypic switching using trophoblast stem cells. They identify Rel, an NF-kB family member, to be transcriptionally coregulated by CBP and HDAC6. Surprisingly, Rel induces mesenchymal-epithelial transition and itself regulated Hdac6 expression and nuclear localization.

Mediator complex component MED13 regulates zygotic genome activation and is required for postimplantation development in the mouse

Understanding factors that regulate zygotic genome activation (ZGA) is critical for determining how cells are reprogrammed to become totipotent or pluripotent. There is limited information regarding how this process occurs physiologically in early mammalian embryos. Here, we identify a mediator complex subunit, MED13, as translated during mouse oocyte maturation and transcribed early from the zygotic genome. Knockdown and conditional knockout approaches demonstrate that MED13 is essential for ZGA in the mouse, in part by regulating expression of the embryo-specific chromatin remodeling complex, esBAF. The role of MED13 in ZGA is mediated in part by interactions with E2F transcription factors. In addition to MED13, its paralog, MED13L, is required for successful preimplantation embryo development. MED13L partially compensates for loss of MED13 function in preimplantation knockout embryos, but postimplantation development is not rescued by MED13L. Our data demonstrate an essential role for MED13 in supporting chromatin reprogramming and directed transcription of essential genes during ZGA.

Characterization of the Fundulus heteroclitus embryo transcriptional response and development of a gene expression-based fingerprint of exposure for the alternative flame retardant, TBPH (bis (2-ethylhexyl)-tetrabromophthalate)

Although alternative Flame Retardant (FR) chemicals are expected to be safer than the legacy FRs they replace, their risks to human health and the environment are often poorly characterized. This study used a small volume, fish embryo system to reveal potential mechanisms of action and diagnostic exposure patterns for TBPH (bis (2-ethylhexyl)-tetrabromophthalate), a component of several widely-used commercial products. Two different concentration of TBPH were applied to sensitive early life stages of an ecologically important test species, Fundulus heteroclitus (Atlantic killifish), with a well-annotated genome. Exposed fish embryos were sampled for transcriptomics or chemical analysis of parent compound and primary metabolite or observed for development and survival through larval stage. Global transcript profiling using RNA-seq was conducted (n = 16 per treatment) to provide a non-targeted and statistically robust approach to characterize TBPH gene expression patterns. Transcriptomic analysis revealed a dose-response in the expression of genes associated with a surprisingly limited number of biological pathways, but included the aryl hydrocarbon receptor signal transduction pathway, which is known to respond to several toxicologically-important chemical classes. A transcriptional fingerprint using Random Forests was developed that was able to perfectly discriminate exposed vs. non-exposed individuals in test sets. These results suggest that TBPH has a relatively low potential for developmental toxicity (at least in fishes), despite concerns related to its structural similarities to endocrine disrupting chemicals and that the early life stage Fundulus system may provide a convenient test system for exposure characterization. More broadly, this study advances the usefulness of a biological testing and analysis system utilizing non-targeted transcriptomics profiling and early developmental endpoints that complements current screening methods to characterize chemicals of ecological and human health concern.

Quantitative analysis of ChIP-seq data uncovers dynamic and sustained H3K4me3 and H3K27me3 modulation in cancer cells under hypoxia

Background

A comprehensive assessment of the epigenetic dynamics in cancer cells is the key to understanding the molecular mechanisms underlying cancer and to improving cancer diagnostics, prognostics and treatment. By combining genome-wide ChIP-seq epigenomics and microarray transcriptomics, we studied the effects of oxygen deprivation and subsequent reoxygenation on histone 3 trimethylation of lysine 4 (H3K4me3) and lysine 27 (H3K27me3) in a breast cancer cell line, serving as a model for abnormal oxygenation in solid tumors. A priori, epigenetic markings and gene expression levels not only are expected to vary greatly between hypoxic and normoxic conditions, but also display a large degree of heterogeneity across the cell population. Where traditionally ChIP-seq data are often treated as dichotomous data, the model and experiment here necessitate a quantitative, data-driven analysis of both datasets.

Results

We first identified genomic regions with sustained epigenetic markings, which provided a sample-specific reference enabling quantitative ChIP-seq data analysis. Sustained H3K27me3 marking was located around centromeres and intergenic regions, while sustained H3K4me3 marking is associated with genes involved in RNA binding, translation and protein transport and localization. Dynamic marking with both H3K4me3 and H3K27me3 (hypoxia-induced bivalency) was found in CpG-rich regions at loci encoding factors that control developmental processes, congruent with observations in embryonic stem cells.

Conclusions

In silico-identified epigenetically sustained and dynamic genomic regions were confirmed through ChIP-PCR in vitro, and obtained results are corroborated by published data and current insights regarding epigenetic regulation.

Electronic supplementary material

The online version of this article (doi:10.1186/s13072-016-0090-4) contains supplementary material, which is available to authorized users.

Expression and methylation data from SLE patient and healthy control blood samples subdivided with respect to ARID3a levels

Previously published studies revealed that variation in expression of the DNA-binding protein ARID3a in B lymphocytes from patients with systemic lupus erythematosus (SLE) correlated with levels of disease activity (“Disease activity in systemic lupus erythematosus correlates with expression of the transcription factor AT-rich-interactive domain 3A” (J.M. Ward, K. Rose, C. Montgomery, I. Adrianto, J.A. James, J.T. Merrill et al., 2014) [1]). The data presented here compare DNA methylation patterns from SLE peripheral blood mononuclear cells obtained from samples with high numbers of ARID3a expressing B cells (ARID3a^H) versus SLE samples with normal numbers of ARID3a⁺ B cells (ARID3a^N). The methylation data is available at the gene expression omnibus (GEO) repository, “Gene Expression Omnibus: NCBI gene expression and hybridization array data repository” (R. Edgar, M. Domrachev, A.E. Lash, 2002) [2]. Isolated B cells from SLE ARID3a^H and ARID3a^N B samples were also evaluated via qRT-PCR for Type I interferon (IFN) signature and pathway gene expression levels by qRT-PCR. Similarly, healthy control B cells and B cells stimulated to express ARID3a with the TLR agonist, CpG, were also compared via qRT-PCR. Primers designed to detect 6 IFNa subtype mRNAs were tested in 4 IFNa, Epstein-Barr Virus-transformed B cell lines (“Reduced interferon-alpha production by Epstein-Barr virus transformed B-lymphoblastoid cell lines and lectin-stimulated lymphocytes in congenital dyserythropoietic anemia type I” (S.H. Wickramasinghe, R. Hasan, J. Smythe, 1997) [3]). The data in this article support the publication, “Human effector B lymphocytes express ARID3a and secrete interferon alpha” (J.M. Ward, M.L. Ratliff, M.G. Dozmorov, G. Wiley, J.M. Guthridge, P.M. Gaffney, J.A. James, C.F. Webb, 2016) [4].

ABSSeq: a new RNA-Seq analysis method based on modelling absolute expression differences

Background

The recent advances in next generation sequencing technology have made the sequencing of RNA (i.e., RNA-Seq) an extemely popular approach for gene expression analysis. Identification of significant differential expression represents a crucial initial step in these analyses, on which most subsequent inferences of biological functions are built. Yet, for identification of these subsequently analysed genes, most studies use an additional minimal threshold of differential expression that is not captured by the applied statistical procedures.

Results

Here we introduce a new analysis approach, ABSSeq, which uses a negative binomal distribution to model absolute expression differences between conditions, taking into account variations across genes and samples as well as magnitude of differences. In comparison to alternative methods, ABSSeq shows higher performance on controling type I error rate and at least a similar ability to correctly identify differentially expressed genes.

Conclusions

ABSSeq specifically considers the overall magnitude of expression differences, which enhances the power in detecting truly differentially expressed genes by reducing false positives at both very low and high expression level. In addition, ABSSeq offers to calculate shrinkage of fold change to facilitate gene ranking and effective outlier detection.

Electronic supplementary material

The online version of this article (doi:10.1186/s12864-016-2848-2) contains supplementary material, which is available to authorized users.

Deficiency of the placenta- and yolk sac-specific tristetraprolin family member ZFP36L3 identifies likely mRNA targets and an unexpected link to placental iron metabolism

The ZFP36L3 protein is a rodent-specific, placenta- and yolk sac-specific member of the tristetraprolin (TTP) family of CCCH tandem zinc finger proteins. These proteins bind to AU-rich elements in target mRNAs, and promote their deadenylation and decay. We addressed the hypotheses that the absence of ZFP36L3 would result in the accumulation of target transcripts in placenta and/or yolk sac, and that some of these would be important for female reproductive physiology and overall fecundity. Mice deficient in ZFP36L3 exhibited decreased neonatal survival rates, but no apparent morphological changes in the placenta or surviving offspring. We found Zfp36l3 to be paternally imprinted, with profound parent-of-origin effects on gene expression. The protein was highly expressed in the syncytiotrophoblast cells of the labyrinth layer of the placenta, and the epithelial cells of the yolk sac. RNA-Seq of placental mRNA from Zfp36l3 knockout (KO) mice revealed many significantly upregulated transcripts, whereas there were few changes in KO yolk sacs. Many of the upregulated placental transcripts exhibited decreased decay rates in differentiated trophoblast stem cells derived from KO blastocysts. Several dozen transcripts were deemed high probability targets of ZFP36L3; these include proteins known to be involved in trophoblast and placenta physiology. Type 1 transferrin receptor mRNA was unexpectedly decreased in KO placentas, despite an increase in its stability in KO stem cells. This receptor is crucial for placental iron uptake, and its decrease was accompanied by decreased iron stores in the KO fetus, suggesting that this intrauterine deficiency might have deleterious consequences in later life.

Identifying differential transcription factor binding in ChIP-seq

ChIP seq is a widely used assay to measure genome-wide protein binding. The decrease in costs associated with sequencing has led to a rise in the number of studies that investigate protein binding across treatment conditions or cell lines. In addition to the identification of binding sites, new studies evaluate the variation in protein binding between conditions. A number of approaches to study differential transcription factor binding have recently been developed. Several of these methods build upon established methods from RNA-seq to quantify differences in read counts. We compare how these new approaches perform on different data sets from the ENCODE project to illustrate the impact of data processing pipelines under different study designs. The performance of normalization methods for differential ChIP-seq depends strongly on the variation in total amount of protein bound between conditions, with total read count outperforming effective library size, or variants thereof, when a large variation in binding was studied. Use of input subtraction to correct for non-specific binding showed a relatively modest impact on the number of differential peaks found and the fold change accuracy to biological validation, however a larger impact might be expected for samples with more extreme copy number variations between them. Still, it did identify a small subset of novel differential regions while excluding some differential peaks in regions with high background signal. These results highlight proper scaling for between-sample data normalization as critical for differential transcription factor binding analysis and suggest bioinformaticians need to know about the variation in level of total protein binding between conditions to select the best analysis method. At the same time, validation using fold-change estimates from qRT-PCR suggests there is still room for further method improvement.

Post-transcriptional regulation of transcript abundance by a conserved member of the tristetraprolin family in Candida albicans

Members of the tristetraprolin (TTP) family of CCCH tandem zinc finger proteins bind to AU-rich regions in target mRNAs, leading to their deadenylation and decay. Family members in Saccharomyces cerevisiae influence iron metabolism, whereas the single protein expressed in Schizosaccharomyces pombe, Zfs1, regulates cell-cell interactions. In the human pathogen Candida albicans, deep sequencing of mutants lacking the orthologous protein, Zfs1, revealed significant increases (> 1.5-fold) in 156 transcripts. Of these, 113 (72%) contained at least one predicted TTP family member binding site in their 3'UTR, compared with only 3 of 56 (5%) down-regulated transcripts. The zfs1Δ/Δ mutant was resistant to 3-amino-1,2,4-triazole, perhaps because of increased expression of the potential target transcript encoded by HIS3. Sequences of the proteins encoded by the putative Zfs1 targets were highly conserved among other species within the fungal CTG clade, while the predicted Zfs1 binding sites in these mRNAs often 'disappeared' with increasing evolutionary distance from the parental species. C. albicans Zfs1 bound to the ideal mammalian TTP binding site with high affinity, and Zfs1 was associated with target transcripts after co-immunoprecipitation. Thus, the biochemical activities of these proteins in fungi are highly conserved, but Zfs1-like proteins may target different transcripts in each species.

The Drosophila Tis11 protein and its effects on mRNA expression in flies

Members of the mammalian tristetraprolin family of CCCH tandem zinc finger proteins can bind to certain AU-rich elements (AREs) in mRNAs, leading to their deadenylation and destabilization. Mammals express three or four members of this family, but Drosophila melanogaster and other insects appear to contain a single gene, Tis11. We found that recombinant Drosophila Tis11 protein could bind to ARE-containing RNA oligonucleotides with low nanomolar affinity. Remarkably, co-expression in mammalian cells with "target" RNAs demonstrated that Tis11 could promote destabilization of ARE-containing mRNAs and that this was partially dependent on a conserved C-terminal sequence resembling the mammalian NOT1 binding domain. Drosophila Tis11 promoted both deadenylation and decay of a target transcript in this heterologous cell system. We used chromosome deletion/duplication and P element insertion to produce two types of Tis11 deficiency in adult flies, both of which were viable and fertile. To address the hypothesis that Tis11 deficiency would lead to the abnormal accumulation of potential target transcripts, we analyzed gene expression in adult flies by deep mRNA sequencing. We identified 69 transcripts from 56 genes that were significantly up-regulated more than 1.5-fold in both types of Tis11-deficient flies. Ten of the up-regulated transcripts encoded probable proteases, but many other functional classes of proteins were represented. Many of the up-regulated transcripts contained potential binding sites for tristetraprolin family member proteins that were conserved in other Drosophila species. Tis11 is thus an ARE-binding, mRNA-destabilizing protein that may play a role in post-transcriptional gene expression in Drosophila and other insects.

Modeling Expression Plasticity of Genes that Differentiate Drug-sensitive from Drug-resistant Cells to Chemotherapeutic Treatment

By measuring gene expression at an unprecedented resolution and throughput, RNA-seq has played a pivotal role in studying biological functions. Its typical application in clinical medicine is to identify the discrepancies of gene expression between two different types of cancer cells, sensitive and resistant to chemotherapeutic treatment, in a hope to predict drug response. Here we modified and used a mechanistic model to identify distinct patterns of gene expression in response of different types of breast cancer cell lines to chemotherapeutic treatment. This model was founded on a mixture likelihood of Poisson-distributed transcript read data, with each mixture component specified by the Skellam function. By estimating and comparing the amount of gene expression in each environment, the model can test how genes alter their expression in response to environment and how different genes interact with each other in the responsive process. Using the modified model, we identified the alternations of gene expression between two cell lines of breast cancer, resistant and sensitive to tamoxifen, which allows us to interpret the expression mechanism of how genes respond to metabolic differences between the two cell types. The model can have a general implication for studying the plastic pattern of gene expression across different environments measured by RNA-seq.

Mobility of the native Bacillus subtilis conjugative plasmid pLS20 is regulated by intercellular signaling

Horizontal gene transfer mediated by plasmid conjugation plays a significant role in the evolution of bacterial species, as well as in the dissemination of antibiotic resistance and pathogenicity determinants. Characterization of their regulation is important for gaining insights into these features. Relatively little is known about how conjugation of Gram-positive plasmids is regulated. We have characterized conjugation of the native Bacillus subtilis plasmid pLS20. Contrary to the enterococcal plasmids, conjugation of pLS20 is not activated by recipient-produced pheromones but by pLS20-encoded proteins that regulate expression of the conjugation genes. We show that conjugation is kept in the default “OFF” state and identified the master repressor responsible for this. Activation of the conjugation genes requires relief of repression, which is mediated by an anti-repressor that belongs to the Rap family of proteins. Using both RNA sequencing methodology and genetic approaches, we have determined the regulatory effects of the repressor and anti-repressor on expression of the pLS20 genes. We also show that the activity of the anti-repressor is in turn regulated by an intercellular signaling peptide. Ultimately, this peptide dictates the timing of conjugation. The implications of this regulatory mechanism and comparison with other mobile systems are discussed.

Bacteria evolve rapidly due to their short generation time and their ability to exchange genetic material, which can occur via different processes, collectively named Horizontal Gene Transfer (HGT). Most bacteria contain, besides a single chromosome, autonomously replicating units called plasmids. Many plasmids carry genes enabling them to be transferred into plasmid-free bacteria. This process, called conjugation, contributes significantly to HGT. Many plasmids also contain antibiotic resistance genes. Therefore, plasmid conjugation plays a major role in the spread of antibiotic resistance. Understanding the regulation of conjugation genes is essential for designing strategies to combat the spread of antibiotic resistance. We have studied the regulation of the native plasmid pLS20 from Bacillus subtilis. Besides being a soil bacterium, B. subtilis is a gut commensal in animals and humans. Here we unraveled the mechanisms controlling conjugation and found that pLS20 conjugation genes become activated when plasmid-free recipient cells are present. We have identified the repressor protein that keeps conjugation in an ‘OFF’ state, and an anti-repressor that activates conjugation. The activity of the anti-repressor is inhibited by a pLS20-encoded peptide that is secreted from the cell and can be absorbed by cells, after a secondary processing step. Ultimately, it is the signaling-peptide that dictates when conjugation genes become activated.

High-resolution Xist binding maps reveal two-step spreading during X-chromosome inactivation

The Xist long noncoding RNA (lncRNA) is essential for X-chromosome inactivation (XCI), the process by which mammals compensate for unequal numbers of sex chromosomes. During XCI, Xist coats the future inactive X chromosome (Xi) and recruits Polycomb repressive complex 2 (PRC2) to the X-inactivation centre (Xic). How Xist spreads silencing on a 150-megabases scale is unclear. Here we generate high-resolution maps of Xist binding on the X chromosome across a developmental time course using CHART-seq. In female cells undergoing XCI de novo, Xist follows a two-step mechanism, initially targeting gene-rich islands before spreading to intervening gene-poor domains. Xist is depleted from genes that escape XCI but may concentrate near escapee boundaries. Xist binding is linearly proportional to PRC2 density and H3 lysine 27 trimethylation (H3K27me3), indicating co-migration of Xist and PRC2. Interestingly, when Xist is acutely stripped off from the Xi in post-XCI cells, Xist recovers quickly within both gene-rich and gene-poor domains on a timescale of hours instead of days, indicating a previously primed Xi chromatin state. We conclude that Xist spreading takes distinct stage-specific forms. During initial establishment, Xist follows a two-step mechanism, but during maintenance, Xist spreads rapidly to both gene-rich and gene-poor regions.

PAVIS: a tool for Peak Annotation and Visualization

We introduce a web-based tool, Peak Annotation and Visualization (PAVIS), for annotating and visualizing ChIP-seq peak data. PAVIS is designed with non-bioinformaticians in mind and presents a straightforward user interface to facilitate biological interpretation of ChIP-seq peak or other genomic enrichment data. PAVIS, through association with annotation, provides relevant genomic context for each peak, such as peak location relative to genomic features including transcription start site, intron, exon or 5'/3'-untranslated region. PAVIS reports the relative enrichment P-values of peaks in these functionally distinct categories, and provides a summary plot of the relative proportion of peaks in each category. PAVIS, unlike many other resources, provides a peak-oriented annotation and visualization system, allowing dynamic visualization of tens to hundreds of loci from one or more ChIP-seq experiments, simultaneously. PAVIS enables rapid, and easy examination and cross-comparison of the genomic context and potential functions of the underlying genomic elements, thus supporting downstream hypothesis generation.

Reconstructing regulatory networks from the dynamic plasticity of gene expression by mutual information

The capacity of an organism to respond to its environment is facilitated by the environmentally induced alteration of gene and protein expression, i.e. expression plasticity. The reconstruction of gene regulatory networks based on expression plasticity can gain not only new insights into the causality of transcriptional and cellular processes but also the complex regulatory mechanisms that underlie biological function and adaptation. We describe an approach for network inference by integrating expression plasticity into Shannon's mutual information. Beyond Pearson correlation, mutual information can capture non-linear dependencies and topology sparseness. The approach measures the network of dependencies of genes expressed in different environments, allowing the environment-induced plasticity of gene dependencies to be tested in unprecedented details. The approach is also able to characterize the extent to which the same genes trigger different amounts of expression in response to environmental changes. We demonstrated the usefulness of this approach through analysing gene expression data from a rabbit vein graft study that includes two distinct blood flow environments. The proposed approach provides a powerful tool for the modelling and analysis of dynamic regulatory networks using gene expression data from distinct environments.

DeconRNASeq: a statistical framework for deconvolution of heterogeneous tissue samples based on mRNA-Seq data

SUMMARY

For heterogeneous tissues, measurements of gene expression through mRNA-Seq data are confounded by relative proportions of cell types involved. In this note, we introduce an efficient pipeline: DeconRNASeq, an R package for deconvolution of heterogeneous tissues based on mRNA-Seq data. It adopts a globally optimized non-negative decomposition algorithm through quadratic programming for estimating the mixing proportions of distinctive tissue types in next-generation sequencing data. We demonstrated the feasibility and validity of DeconRNASeq across a range of mixing levels and sources using mRNA-Seq data mixed in silico at known concentrations. We validated our computational approach for various benchmark data, with high correlation between our predicted cell proportions and the real fractions of tissues. Our study provides a rigorous, quantitative and high-resolution tool as a prerequisite to use mRNA-Seq data. The modularity of package design allows an easy deployment of custom analytical pipelines for data from other high-throughput platforms.

AVAILABILITY

DeconRNASeq is written in R, and is freely available at http://bioconductor.org/packages.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

A quantitative model of transcriptional differentiation driving host-pathogen interactions

Despite our expanding knowledge about the biochemistry of gene regulation involved in host-pathogen interactions, a quantitative understanding of this process at a transcriptional level is still limited. We devise and assess a computational framework that can address this question. This framework is founded on a mixture model-based likelihood, equipped with functionality to cluster genes per dynamic and functional changes of gene expression within an interconnected system composed of the host and pathogen. If genes from the host and pathogen are clustered in the same group due to a similar pattern of dynamic profiles, they are likely to be reciprocally co-evolving. If genes from the two organisms are clustered in different groups, this means that they experience strong host-pathogen interactions. The framework can test the rates of change for individual gene clusters during pathogenic infection and quantify their impacts on host-pathogen interactions. The framework was validated by a pathological study of poplar leaves infected by fungal Marssonina brunnea in which co-evolving and interactive genes that determine poplar-fungus interactions are identified. The new framework should find its wide application to studying host-pathogen interactions for any other interconnected systems.

Genome drafts of four phytoplasma strains of the ribosomal group 16SrIII

By applying a coverage-based read selection and filtration through a healthy plant dataset, and a post-assembly contig selection based on homology and linkage, genome sequence drafts were obtained for four phytoplasma strains belonging to the 16SrIII group (X disease clade), namely Vaccinium Witches' Broom phytoplasma (647 754 nt in 272 contigs), Italian Clover Phyllody phytoplasma strain MA (597 245 nt in 197 contigs), Poinsettia branch-inducing phytoplasma strain JR1 (631 440 nt in 185 contigs) and Milkweed Yellows phytoplasma (583 806 nt in 158 contigs). Despite assignment to different 16SrIII subgroups, the genomes of the four strains were similar, comprising a highly conserved core (92-98 % similar in their nucleotide sequence among each other over alignments about 500 kb in length) and a minor strain-specific component. As far as their protein complement was concerned, they did not differ significantly in their basic metabolism potential from the genomes of other wide-host-range phytoplasmas sequenced previously, but were distinct from strains of other species, as well as among each other, in genes encoding functions conceivably related to interactions with the host, such as membrane trafficking components, proteases, DNA methylases, effectors and several hypothetical proteins of unknown function, some of which are likely secreted through the Sec-dependent secretion system. The four genomes displayed a group of genes encoding hypothetical proteins with high similarity to a central domain of IcmE/DotG, a core component of the type IVB secretion system of Gram-negative Legionella spp. Conversely, genes encoding functional GroES/GroEL chaperones were not detected in any of the four drafts. The results also indicated the significant role of horizontal gene transfer among different 'Candidatus Phytoplasma' species in shaping phytoplasma genomes and promoting their diversity.

GFOLD: a generalized fold change for ranking differentially expressed genes from RNA-seq data

MOTIVATION

RNA-seq has been widely used in transcriptome analysis to effectively measure gene expression levels. Although sequencing costs are rapidly decreasing, almost 70% of all the human RNA-seq samples in the gene expression omnibus do not have biological replicates and more unreplicated RNA-seq data were published than replicated RNA-seq data in 2011. Despite the large amount of single replicate studies, there is currently no satisfactory method for detecting differentially expressed genes when only a single biological replicate is available.

RESULTS

We present the GFOLD (generalized fold change) algorithm to produce biologically meaningful rankings of differentially expressed genes from RNA-seq data. GFOLD assigns reliable statistics for expression changes based on the posterior distribution of log fold change. In this way, GFOLD overcomes the shortcomings of P-value and fold change calculated by existing RNA-seq analysis methods and gives more stable and biological meaningful gene rankings when only a single biological replicate is available.

AVAILABILITY

The open source C/C++ program is available at http://www.tongji.edu.cn/∼zhanglab/GFOLD/index.html

Posttranscriptional regulation of cell-cell interaction protein-encoding transcripts by Zfs1p in Schizosaccharomyces pombe

Members of the tristetraprolin (TTP) family of CCCH tandem zinc finger proteins can bind directly to AU-rich elements in mRNAs and promote transcript deadenylation and decay. The yeast Schizosaccharomyces pombe expresses a single TTP family member, Zfs1p. In this study, we identified probable Zfs1p target mRNAs by comparing transcript levels in wild-type yeast and zfs1Δ mutants, using deep sequencing and microarray approaches. We also used direct RNA sequencing to determine polyadenylation site locations and to confirm the presence of potential Zfs1p target sequences within the target mRNA. These studies identified a set of transcripts containing potential Zfs1p binding sites that accumulated significantly in the zfs1Δ mutants; a subset of these transcripts decayed more slowly in the zfs1Δ mutants and bound directly to Zfs1p in coimmunoprecipitation assays. One apparent direct target encodes the transcription factor Cbf12p, which is known to increase cell-cell adhesion and flocculation when overexpressed. Studies of zfs1Δ cbf12Δ double mutants demonstrated that the increased flocculation seen in zfs1Δ mutants is due, at least in part, to a direct effect on the turnover of cbf12 mRNA. These data suggest that Zfs1p can both directly and indirectly regulate the levels of transcripts involved in cell-cell adhesion in this species.

Differential gene expression in tamoxifen-resistant breast cancer cells revealed by a new analytical model of RNA-Seq data

Resistance to tamoxifen (Tam), a widely used antagonist of the estrogen receptor (ER), is a common obstacle to successful breast cancer treatment. While adjuvant therapy with Tam has been shown to significantly decrease the rate of disease recurrence and mortality, recurrent disease occurs in one third of patients treated with Tam within 5 years of therapy. A better understanding of gene expression alterations associated with Tam resistance will facilitate circumventing this problem. Using a next generation sequencing approach and a new bioinformatics model, we compared the transcriptomes of Tam-sensitive and Tam-resistant breast cancer cells for identification of genes involved in the development of Tam resistance. We identified differential expression of 1215 mRNA and 513 small RNA transcripts clustered into ERα functions, cell cycle regulation, transcription/translation, and mitochondrial dysfunction. The extent of alterations found at multiple levels of gene regulation highlights the ability of the Tam-resistant cells to modulate global gene expression. Alterations of small nucleolar RNA, oxidative phosphorylation, and proliferation processes in Tam-resistant cells present areas for diagnostic and therapeutic tool development for combating resistance to this anti-estrogen agent.