Gene Expression Evolution Figure 2. Relationship between model parameters L, M, N, and exponent b. (A) M is fixed to 300,000. (B) N is fixed to 20,000. Each point shows the exponent b averaged over 50 repetitions of the simulation, where sk = 0.05 and t = 10,000. At the initial state, the mRNA abundance of all genes was set to M/N. The rectangle shows the range of the parameters in typical mammalian cells. doi:10.1371/journal.pone.0007943.g002 affected by M and N. The exponent b takes a value close to 1.0 over a wide range of M and N if and only if L is near 1.0 copy/cell. From this result, the uniformity of exponent b can be explained as follows. The number of expressed genes N and total number of mRNA molecules in a cell M may differ among cells, but those parameters have little effect on the value of b. On the other hand, the lower limit L might be common across a variety of cells, because if fi(t) becomes under 1.0 copy/cell by mutation, that means the transcript can no longer exist at every moment in the cytoplasm. This situation would, in effect, be similar to the gene being lost; hence, it would be subject to negative selection. Therefore L would have a similar value across a variety of cells, and consequently b would be the same across a variety of cells. Loss of Clock-Likeness in the Long-Run Behavior of Gene Expression Divergence under the Neutral Model It is well known that Zipf’s law-like distributions can be generated by a variety of stochastic processes [40,41,46,47]. Those processes include the law of proportionate effect which is described above, the theory of breakage studied by Kolmogoroff [48], the preferential attachment process studied by Yule [49] and Simon [50], the optimization theory of the genesis of Zipf’s law studied by Mandelbrot [51], and the distribution of powers and products of normal variables studied by Haldane [52]. On the other hand, some people even believe the superstition that a Zipf’s law-like distribution can be generated ‘‘without cause’’. Needless to say, this belief is not true. The truth is that in some phenomena the cause of the genesis of the distribution is not known, and in some phenomena the cause is trivial and rather disappointing. Regardless of how many models were introduced previously, the law of proportionate effect is currently a good starting point available for the explanation of the genesis of the Zipf’s law-like distribution in mRNA abundance, since it is the basis for the model of the clock-like accumulation of gene expression divergence [53–55]. It should be recalled that if evolutionary forces affect the mRNA abundance of each gene, they inevitably affect their distributions. Therefore, the explanations of the clock- like accumulation and the Zipf’s law-like distribution must be based on the same set of mechanisms. Both the distribution of mRNA abundance in a cell (A) and the gene expression divergence (B) are solely determined from the gene expression profiles (C), where a gene expression profile is a list of the mRNA abundance of each gene fi(t), where i = 1, …, N, at a given t. Namely, (A) and (B) are two different features of the same source (C). Therefore, if a model insists that it can explain (A) through a formulation of (C), the model should explain (B) simultaneously, and vice versa. One may say that fi(t) can change among species without any effect on the overall distribution of fi(t) in a cell; say f1(0) = 1.0, f2(0) = 2.0 to f1(1) = 2.0, f2(1) = 1.0. Even in such cases, the source of the constraints which force the distribution into being unchanged should be explained from the time-course of fi(t) formulated by the model. Therefore, we examined whether the refined model, which has been introduced to explain the Zipf’s law-like distribution, can also regenerate another independently found phenomenon—the clock- like accumulation of gene expression divergence discovered by Khaitovich et al. (2004) [22]. The Monte Carlo simulation showed that when the previous model (specifically, in the case of L = 0.0, M = 300,000, and N = 20,000) is assumed, it follows that gene expression divergence increases linearly without limit (Figure 3A). On the other hand, the refined model (L = 1.0 copy/cell) predicts that at first gene expression divergence increases linearly, then the increasing rate decreases gradually, and eventually converges to an upper limit near 0.4, although it exhibits substantial stochastic fluctuations (Figure 3B). The upper limit is found to be unaffected by the standard deviation sk of the distribution of kit (Figure 3B), and slightly affected by M and N (data not shown), similar to the case of the genesis of Zipf’s law-like distribution. This prediction of the presence of the upper limit is one of the discriminating features of the refined model which can be used as a criterion for its falsification. The time-course of gene expression divergence reported in Khaitovich et al. (2004) [22] exhibits signs of an upper limit. According to the neutral theory [56] and the generation time effect hypothesis [57–59], the evolutionary rate is a function of the number of generations rather than chronological time. From a human-orangutan comparison, the divergence increases about 0.3 in about 2T = 1.2 million generations, but from human-orangutan and human-macaque, the divergence increases only by about 0.1 PLoS ONE | www.plosone.org 4 November 2009 | Volume 4 | Issue 11 | e7943