Quantitative relationship between synonymous codon usage bias and GC composition across unicellular genomes.

BACKGROUND: Codon usage bias has been widely reported to correlate with GC composition. However, the quantitative relationship between codon usage bias and GC composition across species has not been reported.
RESULTS: Based on an informatics method (SCUO) we developed previously using Shannon informational theory and maximum entropy theory, we investigated the quantitative relationship between codon usage bias and GC composition. The regression based on 70 bacterial and 16 archaeal genomes showed that in bacteria, SCUO = -2.06 * GC3 + 2.05*(GC3)2 + 0.65, r = 0.91, and that in archaea, SCUO = -1.79 * GC3 + 1.85*(GC3)2 + 0.56, r = 0.89. We developed an analytical model to quantify synonymous codon usage bias by GC compositions based on SCUO. The parameters within this model were inferred by inspecting the relationship between codon usage bias and GC composition across 70 bacterial and 16 archaeal genomes. We further simplified this relationship using only GC3. This simple model was supported by computational simulation.
CONCLUSIONS: The synonymous codon usage bias could be simply expressed as 1+ (p/2)log2(p/2) + ((1-p)/2)log2((l-p)/2), where p = GC3. The software we developed for measuring SCUO (codonO) is available at http://digbio.missouri.edu/~wanx/cu/codonO.

Background

All amino acids except Met and Trp are encoded by more than one codon. DNA sequence data from diverse organisms have shown that synonymous codons for any amino acid are not used with equal frequency, even though choices among codons should be equivalent in terms of protein sequences [1-6]. Previous codon usage analyses showed that codon usage bias is very complicated and is associated with various biological factors, such as gene expression level [7-10], gene length [11-13], gene translation initiation signal [14], protein amino acid composition [6,15], protein structure [16,17], tRNA abundance [18-21], mutation frequency and patterns [22,23], and GC composition [24-27]. In this paper, we further explore the relationship between codon usage and GC composition.

GC composition may be described at three levels: 1) Overall GC content. The overall genome GC content in living organisms varies from 25–75% [28]. However, within a single gene, the overall GC content is 7–95%. 2) Local GC composition. Local GC composition is defined based on the positions on the genetic codons. GC1 is the GC composition at the first site of codons, GC2 the GC composition at the second site of codons, and GC3 the GC composition at the third site of codons. 3) The ratio of G/C or A/T within a single strand of DNA. Based on Watson-Crick base pairing rules, the overall GC content is the same between the plus strand and minus strand of the DNA sequence [29]. However, within a single strand, the ratio of G to C and A to T may not be 1. Bacterial genomes were found to be relatively enriched in G over C and T over A, and slightly depleted in G+C, in their weakly evolution-selected positions (intergenic regions and third codon positions) in the leading strand compared with the lagging strand [28]. Although the overall GC content for a genome is reduced by AT-rich intergenic regions of the genome, the gene GC composition is tightly correlated with the genome's overall GC composition. The higher the overall GC content bias, the higher the local GC composition (GC1, GC2, and GC3).

It has been generally accepted that genome GC content is correlated with amino acid usage and codon usage [30]. A very low or very high GC composition is associated with a large codon usage bias. Recently, through a simple mutational model, Knight et al. [31] showed that it was the GC composition that drives codon and amino-acid usage although both mutation and selection play important roles. By using the corresponding analysis of codon usage over 32 bacterial and 8 archaeal genomes, Lynn et al. [32] further showed that codon usage bias was affected by GC composition and environment (e.g., temperature).

However, there has been no systematic, quantitative evaluation of the relationship between codon usage bias and GC composition. In addition, the theoretical basis underlying the relationship between codon usage bias and GC composition has not been illustrated. In this paper, we applied an informatics method [13], which is based on Shannon informatics theory and the entropy theory, to explore the relationship between GC composition and synonymous codon usage bias among 70 bacterial and 16 archaeal genomes. We presented an analytical model to quantify the non-linear relationship between GC3 and a measurement of codon usage bias (synonymous codon usage order, or SCUO), which reveals that GC3 is the key factor driving synonymous codon usage and that this mechanism is independent of species. These results were supported by our simulations. Our results also showed that the asymmetric distribution of G over C and A over T at the third codon position may increase codon usage bias. The underlying mechanisms behind GC composition and codon usage bias are discussed.

Informatics method

We recently developed an informatics method [13] to provide an estimate for the orderliness of synonymous codon usage (SCUO) and the amount of synonymous codon usage bias. This method was based on the Shannon informatics theory and the entropy theory and allows the comparison of codon usage bias within and across genomes.

To calculate SCUO, we created a codon table for the amino acids that have more than one codon, indexed in an arbitrary way, so that we could unambiguously refer to the j-th (degenerate) codon of amino acid i, 1 ≤ i ≤ 18. In mycoplasmas, Trp was also included into the codon table since a standard stop codon TGA encodes Trp in this specific species so that 1 ≤ i ≤ 19. For simplicity, the following description of the method is only based on the standard genetic codon table although the actual SCUO computation considered special cases for different organisms.

Let nrepresent the number of degenerate codons for amino acid i, so 1 ≤ j ≤ n: for example, 1 ≤ j ≤ 6 for leucine, 1 ≤ j ≤ 2 for tyrosine, etc. For each sequence, let xrepresent the number of times that synonymous codon j of amino acid i is present, 1 ≤ i ≤ 18, 1 ≤ j ≤ n. Normalizing the xby their sum over j gives the frequency of the j-th degenerate codon for amino acid i in each sequence.

According to the information theory, we define the entropy Hof the i-th amino acid of the j-th codon in each sequence by

H= -plog p    2

Summing over the codons representing amino acid i gives the entropy of the i-th amino acid in the each sequence

If the synonymous codons for the i-th amino acid were used at random, one would expect a uniform distribution of them as representatives for the i-th amino acid. Thus, the maximum entropy for the i-th amino acid in each sequence is

If only one of the synonymous codons is used for the i-th amino acid, i.e., the usage of the synonymous codons is biased to the extreme, then the i-th amino acid in each sequence has the minimum entropy:

Unlike Shannon's definition of information, Gatlin [33] and Layzer [34] define the information as the difference between the maximum entropy and the actual entropy as an index of orderliness. The greater the information, the more ordered the sequence will be [35]. In our case, this information measures the nonrandomness in synonymous codon usage and therefore describes the degree of orderliness for synonymous codon usage for the i-th amino acid in each sequence.

Let Obe the normalized difference between the maximum entropy and the observed entropy for the i-th amino acid in each sequence, i.e.,

Obviously, 0 ≤ O≤ 1. When the synonymous codon usage for the i-th amino acid is random, O= 0. When this usage is biased to the extreme, 0= 1. Thus, Ocan be thought of as a measure of the bias in synonymous codon usage for the i-th amino acid in each sequence. We designate the statistics Oas the synonymous codon usage order (SCUO) for the i-th amino acid in each sequence.

Let Fbe the composition ratio of the i-th amino acid in each sequence:

Then the average SCUO for each sequence can be represented as

A software package called codonO was written using the C programming language to calculate SCUO for each open reading frame (ORF). This program is available at . The 86 unicellular genomes explored in this paper can be found at our Web site .

Theoretical model between GC3 and codon usage

From the standard genetic code for synonymous codons, each amino acid was encoded by xyU/C and/or xyA/G. The third site of each codon should be purine (R, either U or C) or pyrimidine (Y, either A or G). For those amino acids with double genetic codons, they can only be coded by xyU/C or xyA/G. We define = A, = U, and vice versa. We assume the probability of G+C at the GC3 is p. Then the probability of A+T at the GC3 is 1-p. Let α be P(G), i.e., the probability of G, and β be P(A), i.e., the probability of A. Hence, 0 < α ≤ p, 0 < β ≤ 1 - p. If we apply Equation 2, we will obtain entropy for each amino acid - αlog2α if the third site is G or -(p - α)log2(p - α) if the third position is C (we use base 2 through our theoretical model). The entropy for each amino acid will be - βlog2β if the third site is A or -(1-p-β)log2(1-p-β) if the third position is T. Let us use Gln as an example, which is coded by CGA/G. If the sequence was coded only by this single amino acid, α = p and β = 1-p. Therefore, by applying Equation 3, we can calculate the overall entropy for synonymous codons, -plog2p - (1-p)log2(1-p). In this case, the maximum entropy will be log22 = 1. Thus, the overall codon bias for Gln will be 1 + plog2p + (1-p)log2(1-p). For all of the amino acids with 2 genetic codons, we can use the same equation:

SCUO n2 = 1 + plog2p + (1-p)log2(1-p)    10

We can apply the similar deducing process for amino acids with 3, 4 or 6 genetic codes. Based on our definition of SCUO (Equation 9), the overall codon usage bias will be the sum of each amino acid. Thus, we obtain this first-order approximation for the overall codon bias:

SCUO = 1 + μ*SCUOn2 + ν*SCUOn3 + ω*SCUOn4 + γ*SCUOn6    11

where μ, ν, ω, and γ are associated adjusted weights, which are the amino acid usage frequency of 2, 3, 4 or 6 genetic codes, respectively. The adjusted weights in archaea and bacteria are shown in Table 1.

Table 1

Adjusted factors in archaea and bacteria

	μ	ν	ω	γ
Archaea	0.43	0.08	0.33	0.16
Bacteria	0.42	0.06	0.35	0.17

To simplify the analysis process, we simplify the analytical model as a binary selection model. Among the 18 synonymous codons, only Ile has three genetic codons. All other synonymous codons are encoded in pairs with xyZ and xy . Assuming no additional codon bias results from a biased distribution between A and T or between C and G, i.e., P(A) = P(T), P(G) = P(C), the synonymous codons bias will be simplified to be a binary selection model:

SCUO = 1+ (p/2)log2(p/2) + ((l-p)/2)log2((l-p)/2)     12

In this case, the codon bias is solely determined by GC3, p. This simplified estimation between codon usage bias and GC3, which will reflect the lowest boundary of the codon usage bias for a specific GC3.

Results

GC composition affects SCUO within a single genome

We examined the relationships of SCUO over GC content, GC1, GC2 and GC3 in E. coli (Fig. 1). We only included genes with >200 codons, based on a comparison between gene length and SCUO in E. coli [13] that showed the SCUOs of genes with <200 codons have large fluctuations in codon usage bias. These results demonstrated that the overall content of GC and GC3 in E. coli have the largest impact on SCUO. This was also shown to be universal within other genomes (data not shown). The E. coli genome exhibited three horns (Figs, 1a and 1d). A lower or higher GC or GC3 over the center GC (50.8%) was associated with a relatively higher SCUO. GC1 showed two horns, whereas GC2 did not show this trend.

Figure 1

Relationship between SCUO and GC composition in (a) Relationship between SCUO and the overall GC composition. (b) Relationship between SCUO and GC1. (c) Relationship between SCUO and GC2. (d) Relationship between SCUO and GC3.

Relationship between GC composition and SCUO across archaeal and bacterial genomes

The comparison between GC composition and SCUO across different genomes showed that the non-linear relationships between SCUO and GC3 among the unicellular organisms were independent of species (Figs. 2 and 3). The dependency of SCUO on GC3 was also strong and non-linear. Based on the quadratic regression, in bacteria, SCUO = -2.06 * GC3 + 2.05*(GC3)2 + 0.65, r = 0.91, P < 0.0001. In archaea, SCUO = -1.79 * GC3 + 1.85*(GC3)2 + 0.56, r = 0.89, P < 0.0001. To investigate the impact of GC3 on GC, we combined both archaeal and bacterial data and plotted the relationships between GC3 and GC (Fig. 4a). GC was linearly correlated (GC3 = 1.78*GC-0.35, r = 0.87) with GC3 for 86 unicellular genomes. We also plotted the relationship between GC1 vs. GC3 and GC2 vs. GC3 (Figs. 4b and 4c). We found that GC3 was linearly correlated with both GC1 and GC2. Interestingly, the slope for GC3 vs. GC1 (0.42) was slightly higher than that for GC3 vs. GC2 (0.27).

Figure 2

Correlation between SCUO and GC composition in 16 archaeal genomes. (a) Relationship between SCUO and GC1. (b) Relationship between SCUO and GC2. (c) Relationship between SCUO and GC3, where the red line denotes the binary selection model (Equation 12), the blue line denotes the codon bias resulting from a single amino acid with binary codons (Equation 10) and the green line denotes the mathematical model for 2, 3, 4 and 6 codons (Equation 11), and the maroon line reflects the quadratic regression line SCUO = -1.79 * GC3 + 1.85*(GC3)2 + 0.56, r = 0.89, P < 0.0001. (d) Relationship between SCUO and overall GC.

Figure 3

Correlation between SCUO and GC composition in 70 bacterial genomes. (a) Relationship between SCUO and GC1. (b) Relationship between SCUO and GC2. (c) Relationship between SCUO and GC3, where the red line denotes the binary selection model (Equation 12), the blue line denotes the codon bias resulting from a single amino acid with binary codons (Equation 10) and the green line denotes the mathematical model for 2, 3, 4 and 6 codons (Equation 11), and the maroon line reflects the quadratic regression line SCUO = -2.06 * GC3 + 2.05*(GC3)2 + 0.65, r = 0.91. (d) Relationship between SCUO and overall GC.

Figure 4

The relationship between GC vs. GC3, GC3 vs. GC1, and GC3 vs. GC2. (a) The relationship between GC and GC3 (GC3 = 1.78*GC-0.35, r = 0.87). (b) The relationship between GC3 and GC1 (GC1 = 0.42*GC3 + 0.32, r = 0.82). (c) The relationship between GC3 and GC2 (GC2 = 0.27*GC3 + 0.25, r = 0.73) in unicellular genomes.

We plotted the simplified binary selection model in the GC3 vs. SCUO of both archaea and bacteria. The binary selection model closely fit the lowest boundary of the correlation between GC3 and SCUO. The codon usage bias was generally above the plotted curve because in most cases, P(A) is not equal P(T) and P(G) is not equal to P(C) in the natural world. The distribution of Ile may also increase the codon usage bias.

Simulation of codon usage bias and GC composition

We carried out numerical simulations to further study the relationship between codon usage bias and GC composition. The simulation process was implemented to mimic the unicellular GC composition. For example, the linear relationship of GC3 vs. GC1 (Fig. 5b) and GC3 vs. GC2 (Fig. 5c) in the pseudo ORFs followed those in archaea and bacteria (Figs. 5b and 5c). The simulation showed that the correlation between GC and GC3 (Fig. 5a) was similar to that in archaea and bacteria (Fig. 4a). Moreover, the relationships of GC1 vs. SCUO (Fig. 5d), GC2 vs. SCUO (Fig. 5e), and GC vs. SCUO (Fig. 5g) in the simulation were also similar to those in archaea and bacteria (Figs. 2 and 3). As with the unicellular data, our simulation results demonstrated that there was good agreement between the non-linear relationship of codon usage bias and GC3 (Fig. 5f), and the binary selection model.

Figure 5

Simulation results in the pseudo open reading frames (ORFs). (a) Relationship between the overall GC and GC3. (b) Relationship between GC1 and GC3. (c) Relationship between GC2 and GC3. (d) Relationship between GC and SCUO. (e) Relationship between GC1 and SCUO. (f) Relationship between GC3 and SCUO, where the red line denotes the binary selection model (Equation 12), the blue line denotes the codon bias resulting from a single amino acid with binary codons (Equation 10) and the green line denotes the mathematical model for 2, 3, 4 and 6 codons (Equation 11). (g) Relationship between GC2 and SCUO.

Impact of asymmetric distribution of G/C and A/T on codon usage bias

We defined P(A) = P(T) and P(G) = P(C) for the third codon positions at step 4 in the simulation process. The correlation between codon usage bias and GC3 is shown in Figure 6. Compared with the random distribution of P(A) and P(G), both the SCUO range and SCUO values of the pseduo ORFs are smaller during symmetry distribution of G over C and A over T at the third codon positions. Compared with simulation results shown in Figure 5, this simulation showed that asymmetric distribution of mutational pressure at the third codon positions may increase codon usage bias.

Figure 6

Simulation of symmetry distribution of mutational pressure on codon usage bias. The red line denotes our simplified binary selection model (Equation 12).

Discussion

In this paper, we applied the informatics method [13] to measure SCUO and GC composition within and across unicellular genomes. We first explored the relationship between SCUO and GC composition in E. coli. Our results demonstrated that the relationship between codon usage bias and GC3 formed a "U" shape, which may reflect a direct mathematic consequence of GC3 over codon usage bias. It is interesting that there was a weak horn in the center of "U" shape (Fig. 1), which is similar to previous reports [36]. The genes in the central horns were discussed in Karlin et al. [36]. We can detect some weak horns in some other genomes, such as Lactococcus lactis. However, the central horns were not shown within some other genomes we explored (data not shown). The comparison between GC composition and SCUO confirmed previous reports that GC3 was the most important factor in codon bias among GC, GC1, GC2, and GC3.

The comparison of codon usage and GC composition across 16 archaeal and 70 bacterial genomes showed that a similar, non-linear relationship existed between GC3 and SCUO across the unicellular cells. This also confirmed that GC3 is the main driving force for the codon usage bias. These results were simulated using pseudo ORFs.

To understand the theoretical basis for the relationship between GC composition and codon usage bias, we derived a binary selection model of GC3 based on the Shannon informatics theory and the entropy concept. The non-linear relationship agreed very well with the theoretical equations (Figs. 3 and 4).

Similar to our results with unicellular organisms, the nonlinear relationship between GC3 and codon usage bias has been reported in human and mouse [37]. This indicates that the relationship is independent of species. However, compared to genes in the same species, it is not known whether the features associated with codon usage bias, such as gene expression level, protein structure, etc., are associated with codon usage bias across species. It may be interesting to investigate the relationships between those reported gene features [7-27] and codon usage bias of those horizontally transferred genes in their new adapted hosts. The methodology described here provides a simple way (SCUO) to compare the codon usage bias across species.

GC3 were demonstrated to be tightly associated with cell functional significance [38,39]. Sequence analysis of human receptor tyrosine kinase genes demonstrated that functionally important transmembrane hydrophobic amino acids are specified by codons containing a higher GC frequency at the third bases than are transmembrane neutral amino acids.

We have demonstrated that the asymmetric distribution of G over C and T over A at the third codon positions increased both the codon usage bias values and ranges. The asymmetric distribution of mutational pressure at the third codon positions provides more flexible selection ability during the environmental adaptation process. The asymmetric distribution has been observed in many bacterial genomes, which were found to be relatively enriched in G over C and T over A at GC3 in the leading strand compared with the lagging strand [28]. The effect of asymmetric mutation pressure on the amino-acid composition of proteins has been reported elsewhere [40-42]. Thus, the asymmetric distribution of mutational pressure contributes to the codon usage bias besides the GC3, although we did not include such an effect in our simplified binary selection model.

It is interesting that the GC3 is not asymmetrically distributed (Fig. 2,3). One possible reason might be that the GC content of protein coding region is overall higher than the non-coding region. The other possible reason might be the length of the distribution. Within a single genome, we found that the shorter sequences have a wider distribution of GC3s (data not shown). Within our simulation, we assumed the even distribution of the gene length between 200 and 1000. This might result in the symmetric distribution GC3 in the simulation.

Conclusions

In summary, we developed a simple binary selection model that mimicked the quantitative relationship between the codon usage bias and GC composition in the unicellular organisms, which was supported by systematically characterizing the relationship between codon usage bias and GC composition among 86 unicellular organisms. Our simulation results support this finding by demonstrating that the asymmetric distribution of mutational pressure at the third codon positions has an impact on codon usage bias.

Methods

Genome database

The bacteria and archaea genomic sequences and annotations were downloaded from in August, 2002.

Computational simulation

To simulate the relationship between codon usage bias and GC composition, we generated 15,514 pseudo ORFs with a random size of 200–1000 codons. For each ORF, the nucleotides (A,T,G, or C) for three positions (1s, 2s, and 3s) were generated separately using the following procedure: 1) generate nucleotides for the third codon position randomly with a GC content between 0 and 1; 2) generate nucleotides for the first codon position with a GC content between 0.42*GC3 + 0.18 and 0.42*GC3 + 0.46, which represents a similar linear relationship between GC3 and GC1 in actual genomes, as shown in the Results section; 3) generate nucleotides for the second codon position with a GC content between 0.27*GC3 + 0.11 and 0.27*GC3 + 0.39, which represents a similar linear relationship between GC3 and GC2 in actual genomes; 4) randomly generate A/T and G/C with equal or random probability between A and T, as well as between G and C, for each position based on the associated GC composition; and 5) combine the generated nucleotides as an ORF, and discard any pseudo ORF containing a stop codon. We calculated SCUO for each pseudo ORF using codonO.

Lists of abbreviations used

SCUO: Synonymous Codon Usage Order

GC1: The GC composition for the first site of each codon

GC2: The GC composition for the second site of each codon

GC3: The GC composition for the third site of each codon

Authors' contributions

XFW was responsible for project design and implementation, software implementation, and manuscript writing. DX participated in the development of the binary selection model and co-wrote the manuscript with XFW and JZ. AK and JZ were responsible for development of the informatics method. All authors read and approved the final manuscript.