Automated maximum likelihood separation of signal from baseline in noisy quantal data.

Data recordings often include high-frequency noise and baseline fluctuations that are not generated by the system under investigation, which need to be removed before analyzing the signal for the system's behavior. In the absence of an automated method, experimentalists fall back on manual procedures for removing these fluctuations, which can be laborious and prone to subjective bias. We introduce a maximum likelihood formalism for separating signal from a drifting baseline plus noise, when the signal takes on integer multiples of some value, as in ion channel patch-clamp current traces. Parameters such as the quantal step size (e.g., current passing through a single channel), noise amplitude, and baseline drift rate can all be optimized automatically using the expectation-maximization algorithm, taking the number of open channels (or molecules in the on-state) at each time point as a hidden variable. Our goal here is to reconstruct the signal, not model the (possibly highly complex) underlying system dynamics. Thus, our likelihood function is independent of those dynamics. This may be thought of as restricting to the simplest possible hidden Markov model for the underlying channel current, in which successive measurements of the state of the channel(s) are independent. The resulting method is comparable to an experienced human in terms of results, but much faster. FORTRAN 90, C, R, and JAVA codes that implement the algorithm are available for download from our website.

In vitro recombinants of antibiotic-resistant Chlamydia trachomatis strains have statistically more breakpoints than clinical recombinants for the same sequenced loci and exhibit selection at unexpected loci.

Lateral gene transfer (LGT) is essential for generating between-strain genomic recombinants of Chlamydia trachomatis to facilitate the organism's evolution. Because there is no reliable laboratory-based gene transfer system for C. trachomatis, in vitro generation of recombinants from antibiotic-resistant strains is being used to study LGT. However, selection pressures imposed on in vitro recombinants likely affect statistical properties of recombination relative to naturally occurring clinical recombinants, including prevalence at particular loci. We examined multiple loci for 16 in vitro-derived recombinants of ofloxacin- and rifampin-resistant L(1) and D strains, respectively, grown with both antibiotics, and compared these with the same sequenced loci among 11 clinical recombinants. Breakpoints and recombination frequency were examined using phylogenetics, bioinformatics, and statistics. In vitro and clinical isolates clustered perfectly into two groups, without misclassification, using Ward's minimum variance based on breakpoint data. As expected, gyrA (confers ofloxacin resistance) and rpoB (confers rifampin resistance) had significantly more breakpoints among in vitro recombinants than among clinical recombinants (P < 0.0001 and P = 0.02, respectively, using the Wilcoxon rank sum test). Unexpectedly, trpA also had significantly more breakpoints for in vitro recombinants (P < 0.0001). There was also significant selection at other loci. The strongest bias was for ompA in strain D (P = 3.3 × 10(-8)). Our results indicate that the in vitro model differs statistically from natural recombination events. Additional genomic studies are needed to determine the factors responsible for the observed selection biases at unexpected loci and whether these are important for LGT to inform approaches for genetically manipulating C. trachomatis.

Simplification of reversible Markov chains by removal of states with low equilibrium occupancy.

We present a practical method for simplifying Markov chains on a potentially large state space when detailed balance holds. A simple and transparent technique is introduced to remove states with low equilibrium occupancy. The resulting system has fewer parameters. The resulting effective rates between the remaining nodes give dynamics identical to the original system's except on very fast timescales. This procedure amounts to using separation of timescales to neglect small capacitance nodes in a network of resistors and capacitors. We illustrate the technique by simplifying various reaction networks, including transforming an acyclic four-node network to a three-node cyclic network. For a reaction step in which a ligand binds, the law of mass action implies a forward rate proportional to ligand concentration. The effective rates in the simplified network are found to be rational functions of ligand concentration.

The sequence and analysis of duplication-rich human chromosome 16.

Human chromosome 16 features one of the highest levels of segmentally duplicated sequence among the human autosomes. We report here the 78,884,754 base pairs of finished chromosome 16 sequence, representing over 99.9% of its euchromatin. Manual annotation revealed 880 protein-coding genes confirmed by 1,670 aligned transcripts, 19 transfer RNA genes, 341 pseudogenes and three RNA pseudogenes. These genes include metallothionein, cadherin and iroquois gene families, as well as the disease genes for polycystic kidney disease and acute myelomonocytic leukaemia. Several large-scale structural polymorphisms spanning hundreds of kilobase pairs were identified and result in gene content differences among humans. Whereas the segmental duplications of chromosome 16 are enriched in the relatively gene-poor pericentromere of the p arm, some are involved in recent gene duplication and conversion events that are likely to have had an impact on the evolution of primates and human disease susceptibility.

Predicting phenotype and emerging strains among Chlamydia trachomatis infections.

Chlamydia trachomatis is a global cause of blinding trachoma and sexually transmitted infections (STIs). We used comparative genomics of the family Chlamydiaceae to select conserved housekeeping genes for C. trachomatis multilocus sequencing, characterizing 19 reference and 68 clinical isolates from 6 continental/subcontinental regions. There were 44 sequence types (ST). Identical STs for STI isolates were recovered from different regions, whereas STs for trachoma isolates were restricted by continent. Twenty-nine of 52 alleles had nonuniform distributions of frequencies across regions (p<0.001). Phylogenetic analysis showed 3 disease clusters: invasive lymphogranuloma venereum strains, globally prevalent noninvasive STI strains (ompA genotypes D/Da, E, and F), and nonprevalent STI strains with a trachoma subcluster. Recombinant strains were observed among STI clusters. Single nucleotide polymorphisms (SNPs) were predictive of disease specificity. Multilocus and SNP typing can now be used to detect diverse and emerging C. trachomatis strains for epidemiologic and evolutionary studies of trachoma and STI populations worldwide.

Transformations that preserve detailed balance in Markov models.

Aggregated Markov processes related by similarity transformation are equivalent in that they cannot be distinguished by steady-state experiments. We derive an explicit formula for the set of all detailed-balance preserving similarity transformations between such continuous time Markov chains with N states. The matrices that define the allowed similarity transformations are found to be a simple non-linear function applied to almost any element of the special orthogonal group in N dimensions. Since a model is identifiable only if there is no similarity transformations to an equivalent model, we expect this result to prove useful in the theory of identification of aggregated Markov chains, an enterprise of growing importance as more and more single molecules yield to observation.

Comparative analysis of hepatitis C virus phylogenies from coding and non-coding regions: the 5' untranslated region (UTR) fails to classify subtypes.

BACKGROUND: The duration of treatment for HCV infection is partly indicated by the genotype of the virus. For studies of disease transmission, vaccine design, and surveillance for novel variants, subtype-level classification is also needed. This study used the Shimodaira-Hasegawa test and related statistical techniques to compare phylogenetic trees obtained from coding and non-coding regions of a whole-genome alignment for the reliability of subtyping in different regions. RESULTS: Different regions of the HCV genome yield inconsistent phylogenies, which can lead to erroneous conclusions about classification of a given infection. In particular, the highly conserved 5' untranslated region (UTR) yields phylogenetic trees with topologies that differ from the HCV polyprotein and complete genome phylogenies. Phylogenetic trees from the NS5B gene reliably cluster related subtypes, and yield topologies consistent with those of the whole genome and polyprotein. CONCLUSION: These results extend those from previous studies and indicate that, unlike the NS5B gene, the 5' UTR contains insufficient variation to resolve HCV classifications to the level of viral subtype, and fails to distinguish genotypes reliably. Use of the 5' UTR for clinical tests to characterize HCV infection should be replaced by a subtype-informative test.

Development of forensic assay signatures for ebolaviruses.

Ebolaviruses are a diverse group of RNA viruses comprising five different species, four of which cause fatal hemorrhagic fever in humans. Because of their high infectivity and lethality, ebolaviruses are considered major biothreat agents. Although detection assays exist, no forensic assays are currently available. Here, we report the development of forensic assays that differentiate ebolaviruses. We performed phylogenetic analyses and identified canonical SNPs for all species, major clades and isolates. TaqMan-MGB allelic discrimination assays based on these SNPs were designed, screened against synthetic RNA templates, and validated against ebolavirus genomic RNAs. A total of 45 assays were validated to provide 100% coverage of the species and variants with additional resolution at the isolate level. These assays enabled accurate forensic analysis on 4 "unknown" ebolaviruses. Unknowns were correctly classified to species and variant. A goal of providing resolution below the isolate level was not successful. These high-resolution forensic assays allow rapid and accurate genotyping of ebolaviruses for forensic investigations.

Pattern formation in an N+Q component reaction-diffusion system.

A general N+Q component reaction-diffusion system is analyzed with regard to pattern forming instabilities (Turing bifurcations). The system consists of N mobile species and Q immobile species. The Q immobile species form in response to reactions between the N mobile species and an immobile substrate and allow the Turing instability to occur. These results are valid both for bifurcations from a spatially uniform state and for systems with an externally imposed gradient as in the experimental systems in which Turing patterns have been observed. It is shown that the critical wave number and the location of the instability in parameter space are independent of the substrate concentration. It is also found that the system necessarily undergoes a Hopf bifurcation as the total substrate concentration is decreased. Further, in the case that all the mobile species diffuse at identical rates we show that if the full system is at a point of Turing bifurcation then the N component mobile subsystem is at transition from an unstable focus to an unstable node, and the critical wave number is simply related to the degenerate positive eigenvalue of the mobile subsystem. A sequence of bifurcations that occur in the eigenspectra as the total substrate concentration is decreased to zero is also discussed.

Rapid phylogenetic and functional classification of short genomic fragments with signature peptides.

BACKGROUND: Classification is difficult for shotgun metagenomics data from environments such as soils, where the diversity of sequences is high and where reference sequences from close relatives may not exist. Approaches based on sequence-similarity scores must deal with the confounding effects that inheritance and functional pressures exert on the relation between scores and phylogenetic distance, while approaches based on sequence alignment and tree-building are typically limited to a small fraction of gene families. We describe an approach based on finding one or more exact matches between a read and a precomputed set of peptide 10-mers. RESULTS: At even the largest phylogenetic distances, thousands of 10-mer peptide exact matches can be found between pairs of bacterial genomes. Genes that share one or more peptide 10-mers typically have high reciprocal BLAST scores. Among a set of 403 representative bacterial genomes, some 20 million 10-mer peptides were found to be shared. We assign each of these peptides as a signature of a particular node in a phylogenetic reference tree based on the RNA polymerase genes. We classify the phylogeny of a genomic fragment (e.g., read) at the most specific node on the reference tree that is consistent with the phylogeny of observed signature peptides it contains. Using both synthetic data from four newly-sequenced soil-bacterium genomes and ten real soil metagenomics data sets, we demonstrate a sensitivity and specificity comparable to that of the MEGAN metagenomics analysis package using BLASTX against the NR database. Phylogenetic and functional similarity metrics applied to real metagenomics data indicates a signal-to-noise ratio of approximately 400 for distinguishing among environments. Our method assigns ~6.6 Gbp/hr on a single CPU, compared with 25 kbp/hr for methods based on BLASTX against the NR database. CONCLUSIONS: Classification by exact matching against a precomputed list of signature peptides provides comparable results to existing techniques for reads longer than about 300 bp and does not degrade severely with shorter reads. Orders of magnitude faster than existing methods, the approach is suitable now for inclusion in analysis pipelines and appears to be extensible in several different directions.

The evolutionary rate dynamically tracks changes in HIV-1 epidemics: application of a simple method for optimizing the evolutionary rate in phylogenetic trees with longitudinal data.

Large-sequence datasets provide an opportunity to investigate the dynamics of pathogen epidemics. Thus, a fast method to estimate the evolutionary rate from large and numerous phylogenetic trees becomes necessary. Based on minimizing tip height variances, we optimize the root in a given phylogenetic tree to estimate the most homogenous evolutionary rate between samples from at least two different time points. Simulations showed that the method had no bias in the estimation of evolutionary rates and that it was robust to tree rooting and topological errors. We show that the evolutionary rates of HIV-1 subtype B and C epidemics have changed over time, with the rate of evolution inversely correlated to the rate of virus spread. For subtype B, the evolutionary rate slowed down and tracked the start of the HAART era in 1996. Subtype C in Ethiopia showed an increase in the evolutionary rate when the prevalence increase markedly slowed down in 1995. Thus, we show that the evolutionary rate of HIV-1 on the population level dynamically tracks epidemic events.

Classification of hepatitis C virus and human immunodeficiency virus-1 sequences with the branching index.

Classification of viral sequences should be fast, objective, accurate and reproducible. Most methods that classify sequences use either pair-wise distances or phylogenetic relations, but cannot discern when a sequence is unclassifiable. The branching index (BI) combines distance and phylogeny methods to compute a ratio that quantifies how closely a query sequence clusters with a subtype clade. In the hypothesis-testing framework of statistical inference, the BI is compared with a threshold to test whether sufficient evidence exists for the query sequence to be classified among known sequences. If above the threshold, the null hypothesis of no support for the subtype relation is rejected and the sequence is taken as belonging to the subtype clade with which it clusters on the tree. This study evaluates statistical properties of the BI for subtype classification in hepatitis C virus (HCV) and human immunodeficiency virus-1 (HIV-1). Pairs of BI values with known positive- and negative-test results were computed from 10,000 random fragments of reference alignments. Sampled fragments were of sufficient length to contain phylogenetic signals that grouped reference sequences together properly into subtype clades. For HCV, a threshold BI of 0.71 yields 95.1% agreement with reference subtypes, with equal false-positive and false-negative rates. For HIV-1, a threshold of 0.66 yields 93.5% agreement. Higher thresholds can be used where lower false-positive rates are required. In synthetic recombinants, regions without breakpoints are recognized accurately; regions with breakpoints do not represent any known subtype uniquely. Web-based services for viral subtype classification with the BI are available online.

Evolution of Chlamydia trachomatis diversity occurs by widespread interstrain recombination involving hotspots.

Chlamydia trachomatis is an obligate intracellular bacterium of major public health significance, infecting over one-tenth of the world's population and causing blindness and infertility in millions. Mounting evidence supports recombination as a key source of genetic diversity among free-living bacteria. Previous research shows that intracellular bacteria such as Chlamydiaceae may also undergo recombination but whether this plays a significant evolutionary role has not been determined. Here, we examine multiple loci dispersed throughout the chromosome to determine the extent and significance of recombination among 19 laboratory reference strains and 10 present-day ocular and urogenital clinical isolates using phylogenetic reconstructions, compatibility matrices, and statistically based recombination programs. Recombination is widespread; all clinical isolates are recombinant at multiple loci with no two belonging to the same clonal lineage. Several reference strains show nonconcordant phylogenies across loci; one strain is unambiguously identified as recombinantly derived from other reference strain lineages. Frequent recombination contrasts with a low level of point substitution; novel substitutions relative to reference strains occur less than one per kilobase. Hotspots for recombination are identified downstream from ompA, which encodes the major outer membrane protein. This widespread recombination, unexpected for an intracellular bacterium, explains why strain-typing using one or two genes, such as ompA, does not correlate with clinical phenotypes. Our results do not point to specific events that are responsible for different pathogenicities but, instead, suggest a new approach to dissect the genetic basis for clinical strain pathology with implications for evolution, host cell adaptation, and emergence of new chlamydial diseases.

Genotype turnover by reassortment of replication complex genes from avian influenza A virus.

Reassortment among the RNA segments of Influenza A virus caused the two most recent human influenza pandemics; recently, reassortment has generated viral genotypes associated with outbreaks of avian H5N1 influenza in Asia and Europe. A statistical analysis has been developed for the systematic identification and characterization of reassortant viruses. The analysis was applied to the genes of the replication complex of 152 avian influenza A viruses isolated between 1966 and 2004 from predominantly terrestrial and domestic aquatic avian species. The results indicated that reassortment among these genes was pervasive throughout this period and throughout both the Eurasian and North American lineages of the virus. Evidence is presented that the circulating genotypes of the replication complex are being replaced continually by novel genotypes created by reassortment. No constraints for coordinated reassortment among genes of the replication complex were evident; rather, reassortment almost always proceeded one segment at a time. A maximum-likelihood estimate of the rate of reassortment was derived. For significantly diverged Asian avian influenza A viruses from the period 1991-2004, it was estimated that the median duration between creation of a new genotype and its next segment reassortment was 3 years. Reassortments that introduced previously unobserved influenza genetic material were detected. These findings point to substantial potential for rapid generation of novel avian influenza A viruses, emphasizing the importance of intensive surveillance of these host species in preparation for a possible pandemic.

Polymorphisms in the nine polymorphic membrane proteins of Chlamydia trachomatis across all serovars: evidence for serovar Da recombination and correlation with tissue tropism.

Chlamydia trachomatis is an intracellular bacterium responsible for ocular, respiratory, and sexually transmitted diseases. The genome contains a nine-member polymorphic membrane protein (Pmp) family unique to members of the order Chlamydiales. Genomic and molecular analyses were performed for the entire pmp gene family for the 18 reference serological variants (serovars) and genovariant Ja to identify specific gene and protein regions that differentiate chlamydial disease groups. The mean genetic distance among all serovars varied from 0.1% for pmpA to 7.0% for pmpF. Lymphogranuloma venereum (LGV) serovars were the most closely related for the pmp genes and were also the most divergent, compared to ocular and non-LGV urogenital disease groups. Phylogenetic reconstructions showed that for six of nine pmp genes (not pmpA, pmpD, or pmpE), the serovars clustered based on tissue tropism. The most globally successful serovars, E and F, clustered distantly from the urogenital group for five pmp genes. These pmp genes may confer a biologic advantage that may facilitate infection and transmission for E and F. Surprisingly, serovar Da clustered with the ocular group from pmpE to pmpI, which are located together in the chromosome, providing statistically significant evidence for intergenomic recombination and acquisition of a genetic composition that could hypothetically expand the host cell range of serovar Da. We also identified distinct domains for pmpE, pmpF, and pmpH where substitutions were concentrated and associated with a specific disease group. Thus, our data suggest a possible structural or functional role that may vary among pmp genes in promoting antigenic polymorphisms and/or diverse adhesions-receptors that may be involved in immune evasion and differential tissue tropism.

Using independent open-to-closed transitions to simplify aggregated Markov models of ion channel gating kinetics.

Deducing plausible reaction schemes from single-channel current traces is time-consuming and difficult. The goal is to find the simplest scheme that fits the data, but there are many ways to connect even a small number of states (>2 million schemes with four open and four closed states). Many schemes make identical predictions. An exhaustive search over model space does not address the many equivalent schemes that will result. We have found a canonical form that can express all reaction schemes for binary channels. This form has the minimal number of rate constants for any rank (number of independent open-closed transitions), unlike other canonical forms such as the well established "uncoupled" scheme. Because all of the interconductance transitions in the new form are independent, we refer to it as the manifest interconductance rank (MIR) form. In the case of four open and four closed states, there are four MIR form schemes, corresponding to ranks 1-4. For many models proposed in the literature for specific ion channels, the equivalent MIR form has dramatically fewer links than the uncoupled form. By using the MIR form we prove that all rank 1 topologies with a given number of open and closed states make identical predictions in steady state, thus narrowing the search space for simple models. Moreover, we prove that fitting to canonical form preserves detailed balance. We also propose an efficient hierarchical algorithm for searching for the simplest possible model consistent with a given data set.

The consistent phylogenetic signal in genome trees revealed by reducing the impact of noise.

Phylogenetic trees based on gene repertoires are remarkably similar to the current consensus of life history. Yet it has been argued that shared gene content is unreliable for phylogenetic reconstruction because of convergence in gene content due to horizontal gene transfer and parallel gene loss. Here we test this argument, by filtering out as noise those orthologous groups that have an inconsistent phylogenetic distribution, using two independent methods. The resulting phylogenies do indeed contain small but significant improvements. More importantly, we find that the majority of orthologous groups contain some phylogenetic signal and that the resulting phylogeny is the only detectable signal present in the gene distribution across genomes. Horizontal gene transfer or parallel gene loss does not cause systematic biases in the gene content tree.

Recombination in the genome of Chlamydia trachomatis involving the polymorphic membrane protein C gene relative to ompA and evidence for horizontal gene transfer.

Genome sequencing of Chlamydia trachomatis serovar D has identified polymorphic membrane proteins (Pmp) that are a newly recognized protein family unique to the Chlamydiaceae family. Cumulative data suggest that these diverse proteins are expressed on the cell surface and might be immunologically important. We performed phylogenetic analyses and statistical modeling with 18 reference serovars and 1 genovariant of C. trachomatis to examine the evolutionary characteristics and comparative genetics of PmpC and pmpC, the gene that encodes this protein. We also examined 12 recently isolated ocular and urogenital clinical samples, since reference serovars are laboratory adapted and may not represent strains that are presently responsible for human disease. Phylogenetic reconstructions revealed a clear distinction for disease groups, corresponding to levels of tissue specificity and virulence of the organism. Further, the most prevalent serovars, E, F, and Da, formed a distinct clade. According to the results of comparative genetic analyses, these three genital serovars contained two putative insertion sequence (IS)-like elements with 10- and 15-bp direct repeats, respectively, while all other genital serovars contained one IS-like element. Ocular trachoma serovars also contained both insertions. Previously, no IS-like elements have been identified for Chlamydiaceae. Surprisingly, 7 (58%) of 12 clinical isolates revealed pmpC sequences that were identical to the sequences of other serovars, providing clear evidence for a high rate of whole-gene recombination. Recombination and the differential presence of IS-like elements among distinct disease and prevalence groups may contribute to genome plasticity, which may lead to adaptive changes in tissue tropism and pathogenesis over the course of the organism's evolution.