Published papers

Computer analyses of the entire GenBank database were conducted to examine correlation between splicing sites and codon positions in reading frames. Intron insertion patterns (i.e., splicing site locations with respect to codon positions) have been analyzed for all of the 64 codons of all the eukaryote taxonomic groups: primates, rodents, mammals, vertebrates, invertebrates, and plants. We found that reading frames are interrupted by an intron at a codon boundary (as opposed to the middle of a codon) significantly more often than expected. This observation is consistent with the exon shuffling hypothesis, because exons that end at codon boundaries can be concatenated without causing a frame shift and thus are evolutionarily advantageous. On the other hand, when introns interrupt at the middles of codons, they exist in between the first and second bases much more frequently than between the second and third bases, despite the fact that boundaries between the first and second bases of codons are generally far more important than those between the second and third bases. The reason for this is not clear and yet to be explained. We also show that the length of an exon is a multiple of 3 more frequently than expected. Furthermore, the total length of two consecutive exons is also more frequently a multiple of 3. All the observations above are consistent with results recently published by Long, Rosenberg, and Gilbert (1995).

Free energy values of mRNA tertiary structures around stop codons were systematically calculated to surmise the hairpin-forming potential for all genes in each of the 16 complete prokaryote genomes. Instead of trying to detect each individual hairpin, we averaged the free energy values around the stop codons over the entire genome to predict how extensively the organism relies on hairpin formation in the process of transcription termination. The free energy values of Escherichia coli K-12 shows a sharp drop, as expected, at 30 bp downstream of the stop codon, presumably due to hairpin-forming sequences. Similar drops are observed for Haemophilus influenzae Rd, Bacillus subtilis and Chlamydia trachomatis, suggesting that these organisms also form hairpins at their transcription termination sites. On the other hand, 12 other prokaryotes - Mycoplasma genitalium, Mycoplasma pneumoniae, Synechocystis PCC6803, Helicobacter pylori, Borrelia burgdorferi, Methanococcus jannaschii, Archaeoglobus fulgidus, Methanobacterium thermoautotrophicum, Aquifex aeolicus, Pyrococcus horikoshii, Mycobacterium tuberculosis and Treponema pallidum - show no apparent decrease in free energy value at the corresponding regions. This result suggests that these prokaryotes, or at least some of them, may never form hairpins at their transcription termination sites.

Retrotransposons have generally been known to integrate randomly into host genomes. Jurka, however, showed some consensus sequence patterns at insertion sites of some mammalian retrotransposons and proposed enzymatic involvement that mediates integration. Jurka used about 400 complete human Alu and rodent ID sequences which retain full length with identical flanking sequences at both ends. In our study, more than 25,000 Alu sequences longer than 250 bp were used for comprehensive analysis to identify any consensus sequence(s) preceding the 5' end of Alu elements. 'Entropy' at each nucleotide position within 500 bases of the 5' Alu end was computed. Significant drop of entropy was observed between position -20 and -10 of the 5' end of Alu, suggesting the existence of certain consensus sequence(s) in the region. Frequencies of all possible triplets (total of 64) were measured in the same region. Observation that frequencies of triplets 'aaa,' 'taa,' and 'tta' in the 5' flanking sequences were high is consistent with Jurka. However, frequencies of triplets 'att' and 'aca' were also significantly high, which are not the primary candidates for nick site in Jurka.

Overlapping genes are defined, in this paper, as a pair of adjacent genes whose coding regions are partly overlapping. We systematically analyzed all overlapping genes in the genomes of two closely related species: Mycoplasma genitalium and Mycoplasma pneumoniae. Careful comparisons were made for homologous genes that are overlapped in one species but not in the other. This comparative analysis allows us to propose a model of how overlapping genes emerged in the course of evolution. It was found that overlapping genes were generated primarily due to the loss of a stop codon in either gene, in many cases, the absence of which resulted in elongation of the 3' end of the gene's coding region. More specifically, the loss of the stop codon took place as a result of the following events: deletion of the stop codon (64.4%), point mutation at the stop codon (4.4%), and frame shift at the end of the coding region (6.7%). Overlapping genes, in a sense, can be thought of as the results of evolutionary pressure to minimize genome size. However, our analysis indicates that many overlapping genes, at least in the genomes of M. genitalium and M. pneumoniae, are due to incidental elongation of the coding regions.

Motivation: Genome sequencing projects and further systematic functional analyses of complete gene sets are producing an unprecedented mass of molecular information for a wide range of model organisms. This provides us with a detailed account of the cell with which we may begin to build models for simulating intracellular molecular processes to predict the dynamic behaviour of living cells. Previous work in biochemical and genetic simulation has isolated well-characterized pathways for detailed analysis, but methods for building integrative models of the cell that incorporate gene regulation, metabolism and signaling have not been established. We, therefore, were motivated to develop a software environment for building such integrative models based on gene sets, and running simulations to conduct experiments in silico. Results: E-CELL, a modeling and simulation environment for biochemical and genetic processes, has been developed. The E-CELL system allows a user to define functions of proteins, protein-protein interactions, protein-DNA interactions, regulation of gene expression and other features of cellular metabolism, as a set of reaction rules. E-CELL simulates cell behavior by numerically integrating the differential equations described implicitly in these reaction rules. The user can observe, through a computer display, dynamic changes in concentrations of proteins, protein complexes and other chemical compounds in the cell. Using this software, we constructed a model of a hypothetical cell with only 127 genes sufficient for transcription, translation, energy production and phospholipid synthesis. Most of the genes are taken from Mycoplasma genitalium, the organism having the smallest known chromosome, whose complete 580 kb genome sequence was determined at TIGR in 1995. We discuss future applications of the E-CELL system with special respect to genome engineering. Availability: The E-CELL software is available upon request. Supplementary information: The complete list of rules of the developed cell model with kinetic parameters can be obtained via our web site at: http://e-cell.org/. Contact: mt@@@sfc.keio.ac.jp.

The translation initiation mechanism of archaebacteria is still not clearly understood. Our previous work showed that ATG triplets before start codons have been strongly depleted in eukaryotic genomes, presumably because ribosome of eukaryotes scans mRNA from the 5' to 3' direction to find proper start codons. Extra ATG triplets before start codons would confuse the process and thus they have been negatively selected in eukaryotic genomes. In eubacterial genomes, on the other hand, ribosome binds to the Shine-Dalgarno (SD) sequence at once without mRNA scanning, and the characteristic patterns of ATG triplet depletion were not observed (Saito, R., Tomita, M., 1999. On negative selection against ATG triplets near start codons in eukaryotic and procaryotic genomes. J. Mol. Evol. 48, 213-217). The ATG triplet analysis on archaebacterial genomes revealed that Methanococcus jannaschii and Pyrococcus horikoshii show patterns similar to eukaryotes, implying that these species employ scanning of mRNA from the 5' to 3' direction in the process of translation initiation. On the other hand, our earlier study found that these archaea have SD-like sequences, which are complementary to the 3' end sequence of 16S rRNA, as in eubacterial translation initiation (Osada, Y., Saito, R., Tomita, M. Analysis of base-pairing potentials between 16S rRNA and 5' UTR for translation initiation in various procaryotes. Bioinformatics, in press). These two results combined lead us to conclude that these archaea probably use a hybrid mechanism; their ribosome scans mRNAs from the 5' to 3' direction and then 16S rRNA binds to the SD-like sequence of the 5' UTR. (C) 1999 Elsevier Science B.V. All rights reserved.

Motivation: It is well accepted that the 3' end of 16S rRNA is directly involved in prokaryotic translation initiation by pairing with the Shine-Dalgarno (SD) sequence, which is located in the ribosome-binding site of mRNA. According to Shine and Dalgarno, Escherichia coli's 5' UTR has the pattern of 'AGGAGG' (SD sequence), which is complementary to the 3' end sequence of 16S rRNA. In this work, we systematically calculated free-energy values of the base pairing between the 3' end of 16S rRNA and the 5' UTR of mRNA, in order to analyze the base-pairing potentials in various prokaryotes. The free-energy values were then plotted over distances from the start codon to visualize the free-energy pattern of 5' UTRs. Results: The average free-energy values fell sharply before the start codon in E. coli, which is consistent with the model that the 3' end of 16S rRNA base pairs with the SD sequence. Haemophilus influenzae, Bacillus subtilis and Helicobacter pylori show a similar pattern, suggesting that the organisms have basically the same mechanism of translation initiation as E. coli. Other eubacteria, such as Synechocystis PCC6803, Mycoplasma genitalium, Mycoplasma pneumoniae and Borrelia burgdorferi also show decreases in their free-energy values, although they are less evident. We also did the same analysis with a eukaryote genome as a control; no fall in free-energy values was observed between the 3' end of 18S rRNA and 5' UTRs of Saccharomyces cerevisiae, suggesting that this organism does not base pair in translation initiation. The three archaebacteria A. fulgidus, M. jannaschii and M. thermoautotrophicum show patterns similar to eubacteria, but not to S. cerevisiae, indicating that archaebacteria are closer to eubacteria than to eukaryotes with respect to the mechanism of translation initiation. From these observations, it appears that the shape of the curve produced by the algorithm can be used to predict the mechanism of translation initiation. Availability: The C programs used in our analysis are available upon request.

Computer analyses of various genome sequences revealed the existence of certain periodical patterns of adenine-adenine dinucleotides (ApA). For each genome sequence of 13 eubacteria, 3 archaebacteria, 10 eukaryotes, 60 mitochondria, and 9 chloroplasts, we counted frequencies of ApA dinucleotides at each downstream position within 50 bp from every ApA. We found that the complete genomes of all three archaebacteria have clear ApA periodicities of about 10 bps. On the other hand, all of the 13 eubacteria we analyzed were found to have an ApA periodicity of about 11 bp. Similar periodicities exist in the 10 eukaryotes, although higher organisms such as primates tend to have weaker periodic patterns. None of the mitochondria and chroloplasts we analyzed showed an evident periodic pattern.

In this study, we analyzed the correlation between codon usage bias and Shine-Dalgarno (SD) sequence conservation, using complete genome sequences of nine prokaryotes. For codon usage bias, we adopted the codon adaptation index (CAI), which is based on the codon usage preference of genes encoding ribosomal proteins, elongation factors, heat shock proteins, outer membrane proteins, and RNA polymerase subunit proteins. To compute SD sequence conservation, we used SD motif sequences predicted by Tompa and systematically aligned them with 5′UTR sequences. We found that there exists a clear correlation between the CAI values and SD sequence conservation in the genomes of Escherichia coli, Bacillus subtilis, Haemophilus influenzae, Archaeoglobus fulgidus, Methanobacterium thermoautotrophicum, and Methanococcus jannaschii, and no relationship is found in M. genitalium, M. pneumoniae, and Synechocystis. That is, genes with higher CAI values tend to have more conserved SD sequences than do genes with lower CAI values in these organisms. Some organisms, such as M. thermoautotrophicum, do not clearly show the correlation. The biological significance of these results is discussed in the context of the translation initiation process and translation efficiency.

Study of the cell will never be complete unless its dynamic behavior is understood. The complex behavior of the cell cannot be determined or predicted unless a computer model of the cell is constructed and computer simulation is undertaken. Rapid accumulation of biological data from genome, proteome, transcriptome and metabolome projects can bring us to the point where it is no longer purely speculative to discuss how to construct virtual cells in silico. This article describes attempts to construct whole cell models. The E-CELL project has completed a couple of virtual cell models, and computer simulations have revealed some biological surprises.

The RIKEN Mouse Gene Encyclopaedia Project, a systematic approach to determining the full coding potential of the mouse genome, involves collection and sequencing of full-length complementary DNAs and physical mapping of the corresponding genes to the mouse genome. We organized an international functional annotation meeting (FANTOM) to annotate the first 21,076 cDNAs to be analysed in this project. Here we describe the first RIKEN clone collection, which is one of the largest described for any organism. Analysis of these cDNAs extends known gene families and identifies new ones.

We conducted comprehensive sequence analysis of 5' flanking regions of primate Alu elements. Information contents were computed and frequencies of 1024 pentanucleotides were measured to approximate the location of a characteristic sequence and to specify its pattern(s), which may be involved in the integration of Alu elements into their host genomes. A large number of samples was used, the wide region of the 5' end of Alu elements was analyzed, and comparisons were made among different subfamilies. Through our analyses, 'TTTTAAAAA' or '(T)(ms)(A)(n)' can be stated as a candidate for the characteristic sequence pattern, which resides around the region 5 to 20 base pairs upstream of the 5' end of Alu elements. This characteristic sequence pattern was more prominent in the sequences of younger Alus, which is a strong indication that the sequence pattern has a role at the time of Alu integration.

The Chi sequences are specific oligomers that stimulate DNA repair by homologous recombination, and are different sequences in each organism. Approximately 75% of the copies of the Chi sequence (5′-GCTGGTGG-3′) of Escherichia coli reside on the leading strand, and this orientation bias is often believed to be a consequence of the biological role of Chi sequences as the signal sequence of RecBCD pathway in DNA replication. However, our computer analysis found that many G-rich oligomers also show this asymmetric orientation pattern. The shift in the Chi orientation bias appears around the replication origin and terminus, but these locations are also coincident with the shift points in GC content or GC skew. We conducted the same analysis with the genome of Bacillus subtilis, and found that in addition to Chi, other G-rich oligomers show similar asymmetric orientation patterns, whose shift points were coincident with those of the GC skew. However, the genome of Haemophilus influenzae Rd, whose GC skew is not so pronounced, does not clearly show asymmetric orientation patterns of Chi or other G-rich oligomers. These results lead us to suggest that the uneven distribution of the Chi orientation between the two strands of the double helix is mostly due to the uneven distribution of G content (GC skew) and that the replication-related function of Chi sequences is not the primary factor responsible for the evolutionary pressure causing the orientation bias. © 2000 Elsevier Science B.V.

The codon adaptation index (CAI) values of all protein-coding sequences of the full-length cDNA libraries of Mus musculus were computed based on the RIKEN mouse full-length cDNA library. We have also computed the extent of consensus in flanking sequences of the initiator ATG codon based on the 'relative entropy' values of respective nucleotide positions (from -20 to +12 bp relative to the initiator ATG codon) for each group of genes classified by CAI values. With regard to the two nucleotides positions (-3 and +4) known to be highly conserved in Kozak's consensus sequence, a clear correlation between CAI values and relative entropy values was observed at position -3 but this was not significant at position +4, although a significant correlation was found at position -1 of the consensus sequence. Further, although no correlation was observed at any additional positions, relative entropy values were very high at positions -4, -6, and -8 in genes with high CAI values. These findings suggest that the extent of conservation in the flanking sequence of the initiator ATG codon including Kozak's consensus sequence was an important factor in modulation of the translation efficiency as well as synonymous codon usage bias particularly in highly expressed genes. © 2001 Published by Elsevier Science B.V.

Background: The HOX gene family consists of highly conserved transcription factors that specify the identity of the body segments along the anteroposterior axis of the embryo. Because the phenotypes of mice with targeted disruptions of Hox genes resemble some patterns of human malformations, mutations in HOX genes have been expected to be associated with a significant number of human malformations. Thus far, however, mutations have been documented in only three of the 39 human HOX genes (HOXD13, HOXA13, and HOXA11) partly because current knowledge on the complete coding sequence and genome structure is limited to only 20 of the 39 human HOX genes. Methods: Taking advantage of the human and mouse draft genome sequences, we attempted to characterize the remaining 19 human HOX genes by bioinformatic analysis including phylogenetic footprinting, the probabilistic prediction method, and comparison of genomic sequences with the complete set of the human anonymous cDNA sequences. Results: We were able to determine the full coding sequences of 19 HOX genes and their genome structure and successfully designed a complete set of PCR primers to amplify the entire coding region of each of the 39 HOX genes from genomic DNA. Conclusions: Our results indicate the usefulness of bioinformatic analysis of the draft genome sequences for clinically oriented research projects. It is hoped that the mutation panel provided here will serve as a launchpad for a new discourse on the genetic basis of human malformations. © 2002 Wiley-Liss, Inc.

Although many studies on alternative splicing of specific genes have been reported in the literature, the general mechanism that regulates alternative splicing has not been clearly understood. In this study, we systematically aligned each pair of the 21,076 cDNA sequences of Mus musculus, searched for putative alternative splicing patterns, and constructed a list of potential alternative splicing sites. Two cDNAs are suspected to be alternatively spliced and originating from a common gene if they share most of their region with a high degree of sequence homology, but parts of the sequences are very distinctive or deleted in either cDNA. The list contains the following information: (1) tissue, (2) developmental stage, (3) sequences around splice sites, (4) the length of each gapped region, and (5) other comments. The list is available at http://www.bioinfo.sfc.keio.ac.jp/ intron. Our results have predicted a number of unreported alternatively spliced genes, some of which are expressed only in a specific tissue or at a specific developmental stage. ©2002 by Cold Spring Harbor Laboratory Press.

We conducted comprehensive analyses on intron positions in the Mus musculus genome by comparing genomic sequences in the GenBank database and cDNA sequences in the mouse cDNA library recently developed by Riken Genomic Sciences Center. Our results confirm that introns have a tendency to be located toward the 5′ end of the gene. The same type of analysis was conducted in the coding region of seven eukaryotes (Saccharomyces cerevisiae, Plasmodium falciparum, Caenorhabditis elegans, Drosophila melanogaster, M. musculus, Homo sapiens, Arabidopsis thaliana). Introns in genes with a single intron have a locational bias toward the 5′ end in all species except A. thaliana. We also measured the distance from the start codon to the position of the intron, and found that single introns prefer the location immediately after the start codon in S. cerevisiae and P. falciparum. We discuss three possible explanations for these findings: (1) they are the consequence of intron loss by reverse-transcriptase; (2) they are necessary to accommodate the function; and (3) they are concerned with the mechanism of pre-mRNA splicing. © 2002 Elsevier Science B.V. All rights reserved.

We describe a method, based on pressure-assisted capillary electrophoresis coupled to electrospray ionization mass spectrometry (PACE/ESI-MS), that allows the simultaneous and quantitative analysis of multivalent anions, such as citrate isomers, nucleotides, nicotinamideadenine dinucleotides, and flavin adenine dinucleotide, and coenzyme A (CoA) compounds. Key to the analysis was using a noncharged polymer, poly(dimethylsiloxane), coated to the inner surface of the capillary to prevent anionic species from adsorbing onto the capillary wall. It was also necessary to drive a constant liquid flow toward the MS by applying air pressure to the inlet capillary during electrophoresis to maintain a conductive liquid junction between the capillary and the electrospray needle. Although theoretical plates were inferior to those obtained by CE/ESI-MS using a cationic polymer-coated capillary, the PACE/ESI-MS method improved reproducibility and sensitivity of these anions. Eighteen anions were separated by PACE and selectively detected by a quadrupole mass spectrometer with a sheath-flow electrospray ionization interface. The relative standard deviations (n = 6) of the method were better than 0.6% for migration times and between 1.4% and 6.2% for peak areas. The detection limits for these species were between 0.4 and 3.7 μmol/L with pressure injection of 50 mbar for 30 s (30 nL), that is, mass detection limits calculated in the range from 12 to 110 fmol at a signal-to-noise ratio of 3. The utility of the method was demonstrated by analysis of Citrate isomers, nucleotides, dinucleotides, and CoA compounds extracted from Bacillus subtilis cells.

Connectivity between 30 cytoarchitectural areas covering almost the whole region of neocortex has been proposed. The data source provided by J.L. Conel was used for the analyses. The analyses consisted of computing correlation coefficients between each cortical layer of one cortical area to each cortical layer of another cortical area over the eight age-points and comparison of the time of the first appearance of myelinated axons per layer per area. Derived connectivity supports some interesting findings reported recently by others and is in good agreement with neural network model deduced from neuroimaging and psychological experiments. (C) 2000 Elsevier Science B.V. All rights reserved.