pnas-1955-kimura-solution of a process of random genetic drift with a continuous model

GENETICS: MOTOO KIMURA 28See, for example, C. W. Beck, R. G. Stevenson, Jr., and L. LaPaz, Popular Astronomy, 59, No. 2 (February, 1951), 93 if. 29 W. Wahl, Geochim. et cosmochim. acta, 2, 91, 1952. 30 See H. C. Urey, Geochim. et cosmochim. acta, 1, 209-277, 1951; see, particularly, pp. 233 ff. This paper was reproduced in The Planets, pp. 124 and 141, particularly. Later G. P. Kuiper published the same idea in Atmospheres of the Earth and Planets (Chicago: University of Chicago Press, 1952), pp. 330, 314 (ref. 18), and 421. 31 See H. Jeffreys, The Earth (4th ed.; Cambridge: At the University Press, 1952), chap iv, p. 159, and papers there cited; H. C. Urey, The Planets, pp. 20 ff. 32 See H. C. Urey, The Planets, p. 124, and Astrophys. J., Suppl., 1, 155, 1954. 33 H. C. Urey, The Planets, pp. 49 ff. 34 See B. Gutenberg, Internal Constitution of the Earth (New York: Dover Publications, Inc., 1951), pp. 150 ff. 36 H. C. Urey, Proc. Roy. Soc. London, A, 219, 281,1953. 36 H. C. Urey, The Planets, p. 176. 37 See H. C. Urey, Phys. Rev., 80, 295,1950; The Planets, pp. 64 ff. n D. B. McLaughlin, Observatory, 74, 166, 1954; Pub. Astron. Soc. Pac., 66, 161, 221, 1954. 39 M. Davidson, J. Brit. Astron. Assoc., 29, 194, 1919; see H. C. Urey, The Planets, pp. 43 ff. 4 H. C. Urey, Astrophys. J., Suppl., 1, 167-172, 1954. 41 A. Blaauw, Bull. Astron. Inst. Neth., 11, 405, 414, 1952; A. Blaauw and M. P. Savendoff, Bull. Astron. Inst. Neth., 12, 69, 1953. 42 J. R. Simonton, R. A. Rightmire, A. L. Long, and T. P. Kohman, unclassified report No. NYO-66.25; Phys. Rev. (in press). SOLUTION OF A PROCESS OF RANDOM GENETIC DRIFT WITH A CONTINUOUS MODEL* BY MOTOO KIMURA DEPARTMENT OF GENETICS, t UNIVERSITY OF WISCONSIN Communicated by Sewall Wright, November 15, 1954 The problem of random genetic drift in finite populations due to random sampling of gametes in reproduction was first treated mathematically by R. A. Fisher,1 using a differential equation. Fisher's general method was appropriate, but owing to the omission of a term in the equation, his result for the rate of decay of variance was only half large enough. The correct solution for the state of steady decay was first supplied by S. Wright,2 using the method of path coefficients and an integral equation. Later Fisher3 corrected his results and also elaborated the terminal part of the distribution in the statistical equilibrium by his method of functional equations. In all these works, however, it was assumed that a state of steady decay had been attained, but nothing was known about the complete solution which might show how the process finally leads to the state of steady decay. The present writer,4 by cal- culating the moments of the distribution and with the help of the Fokker-Planck equation, obtained a solution which assumed an infinite series under the continuous model, showing that the process approaches asymptotically the state of steady decay. At that time, however, only the first few coefficients in the terms of the series could be determined. Pursuing the pioblem further, he arrived at the com- plete solution, which will be reported here. After obtaining these results, the writel recently discovered the work of S. Goldberg.5 In his unpublished thesis Goldberg 144 PROC. N. A. S. GENETICS: MOTOO KIMURA I .. solved the diffusion equation for the gene-frequency distribution in a finite popu- lation when recurrent mutations occur. His solutions have a direct connection with the frequency distribution of unfixed classes in the case of pure random genetic drift, as will be seen below. Consider a random mating population of N diploid parents. Let A and A' be a pair of alleles with frequencies x and 1 - x, respectively. In order to single out the effect of random drift, we shall assume an idealized situation in which selection, migration, and mutation are absent and generations do not overlap. The process of the change in gene frequencies is most adequately described by giving the fre- quency distribution f(x, t) at the tth generation, where x takes on a series of dis- crete values: 0, 1/2N, 2/2N, .. . , 1 - 1/2N, 1. For fairly large N, however, x can be treated as a continuous variable without serious error. First, we shall derive the moment formula which is useful heuristically. Let xi be the gene frequency in the tth generation, and let 8t, be the amount of change due to random sampling of gametes per generation, such that XI+1 = Xi + ax,. (1) Let us('+1) = E(xt+) be the nth moment of the distribution about 0 in the (t + 1)th generation. We write the expectation of Xn4i in terms of (x: + 5z1)n in two steps: (1) taking expectation for the random change a5x, which will be de- noted by Ea, and (2) taking the expectation for the existing distribution, which will be denoted by E+. Noting that EB(6x,) = 0, E(6x,))2 = XI(1 - x,)/2N, etc., nl,+1) = E(Xg + 8X,)' E,{xOn + ()Xn-'1Es (ax) + ()x -2Es(bxt)2 + * * (2) E { s+ n(n- 1) z-2 XI(2N X) 0 The intrinsic assumption in the continuous model is that the effective size N is sufficiently large so that terms of -order 1/N2 and higher can- be omitted without serious error. Thus 5)={1 -)}+n(n 1)n--*) For large N, the moments change very slowly per generation, and we -can replace the above equation by the system of differential equations: __ -1=_ n(8 , I) A I I (n 1, 2,3,.. ). (3) If the population starts from the.gene frequency p (0 < p < 1), , ni$°)= pfl and we can obtain the nth moment as the solution of (3), co ,A = p + (2i + 1)pq(-1)F(1 -i, i +22, p) X (n +1)... (n + )(n-1)*- (n-i)4¢[i(i+114f 4(n + 1)...(n + i)v VOL. 41, 1955 145 GENETICS: MOTOO KIMURA where F(1- i, i + 2, 2, p) is the hypergeometric function and q = 1- p. For finite n the series is finite. Putting n = 1, 2, 3, and 4, it will be seen that, for large N, the resulting formulas give a very good approximation to the exact moment formulas obtained by A. Robertson6 (p. 205). The probability f(1, t) of the gene A being fixed in the population by the tth generation can be obtained by using the relation .~~~~~~~~~~~ f(1, t) = lim E xnf(x, t) = lim t, n - co x = n co The resulting series, f(1, t) = p + E (2i + 1)pq(-1)F(1 - i, i + 2, 2, p)e-[t(i+l)/4N]t, (5) is now an infinite series whose convergence must be examined. It is convenient here to introduce the Gegenbauer polynomial Tj!1(z) which is related to the hyper- geometric function by i(i+) F i21-, 1 Z)TiIW 2 ( Y1 2 2) The properties of this function have been thoroughly studied (see Morse and Feshbach7 [pp. 782-783]). Using this relation and putting p = (1 - r)/2, (-1 < r < 1), we obtain "O (2i+1)114if(1, t) = P + E 1 2i(i + 1) (1 - r2) TI-1 (r)e&[i(i+)/NlI (6) where Tol(r) = 1, T11(r) = 3r, T21(r) = 3/2(5r2 - 1), T3I(r) = 5/2(7r3- 3r), etc. Here, if we use the recurrence relation: (2i + 1)(1 - r2)T.li(r) = i(i + 1)P1_(r) - i(i + 1)Pi+i(r), the above formula becomes f1 ) = p EI_+{p i(r) - Pj+,(r) I + ) (7) i= 2 where Pn(r) represents a Legendre polynomial: Po = 1, P1 = r, P2 = 1/2(3r2 - 1), P3 = 1/2(5r2- 3r), etc. For t = 0, the partial sum of the first n terms of equation (7) is (- 1)n-i(Pn1 - Pn)/2 (n > 3). By using a proper integral expression (see later part), it can easily be shown that, if -1 < r < 1, the partial sum tends to zero as n goes to infinity. For t > 0, the series (7) is uniformly convergent and tends to p as t -- c. The probability of the gene A' being fixed in the population by the tth generation f(0, t) can be obtained by replacing p with q and r with -r. If we note that Pn((-r) = (-1)nPP(r), the frequency of the fixed classes is seen to be f(1 t) + f(0, t) = 1I-{P2J(r) - P2J+2(r) } e-[(2j+l)(2j+2)/4NJt (8) which is 0when t = 0 and tends to 1 as t *. Let us now consider the probability distribution of unfixed classes. Let 4(x, t) be the probability density that the gene frequency in the tth generation is between x and x + dx (0 < x < 1). It has been shown that, under the assumption of the continuous model, qb(x, t) satisfies the following partial differential equation: 146 PROC. N. A. S. GENETICS: MOTOO KIMURA 1 a2 x)0, (9) which is a Fokker-Planck equation for the case of the random drift.4' 8 This equation has singularities at the boundaries, and no arbitrary conditions can be imposed there.9 But the moment formula which can be obtained from equations (4) and (5) by calculating, -M_ n'f(1, t) = folxno(x, t) dX, (10) suggests that equation (9) must have the solution of the form co E CjXj(x)e [i(i + 1)/4N]t wheie C1 are constants and X1(x) are functions of x only. In order to solve equa- tion (9), if we put 4 a: Xi(x) exp {-i(i + 1)t/4N} (i = 1, 2, 3, ...), we obtain the hypergeometric equation d2X1 dX1 x(1-x) dX + 2(1-2x)--(1-i)(i + 2)Xi0=OdX2 dx or, putting x = (1 - z)/2 such that z = 1 - 2x (-1 < z < 1), we obtain the Gegen- bauer equation: d2X dXj(z2 - 1) d-2 + 4zd-- (i - 1)(i + 2)X1 = 0. (11) From the comparison of the results obtained from equation (10), it can be found that a Gegenbauer polynomial X1 = T>_(z) is a pertinent solution. Thus cD +(x, t) = E CiTi._i(z)e-i(i + 1)/4Nf]t (12) i= 1 The coefficients C( can be determined from the initial condition that the popula- tion starts from the gene frequency p. Mathematically, (x -p) = E C1T1!(z), (13) i = 1 where S(x) represents the delta function. Multiplying (1 - z2)T'l-(z) on both sides of equation (13) and using the orthogonal property, 2(i + 1)i -1(z2)T (Z)T 1(z) dI = 2m, t-1 (2i + 1) (14) where m in Kronecker's notation.represents zero or a positive integer, we obtain 2{1 - (1 - 2p)2}Tl_(I - 2p) C 2(i + 1)i(2i +1) or (2i + 1) 1Ci =4pq Tj i(1- 2p).i~+) VOL. 41, 1955 147 GENETICS: MOTOO KIMURA Thus the formal solution is co (2i + 1) (1-r) 1__- [ (15) i= 1 i(i+ 1) (5 or, in terms of the hypergeometric function, co O(x, t) = Z pqi(i + 1) (2i + 1)F(1- ii + 2, 2, p) X F(1- i, i + 2, 2, x)e-i(i+l)/4Nlt (15)' For t > 0, the series is uniformly convergent for x and p, since the exponential term approachs zero rapidly. It is interesting that this solution agrees with the "absorb- ing barrier solution" of Goldberg5 by putting a = P = 0 in his formula. The probability that both A and A' coexist in the population in the tth generation (%) is easily obtained from equation (15), by noting that dP,(z)/dz = Tl-(z) and PnM() = 1: := f 4(x, t) dx= (x t) dz o -1 ~~~~2 = (4m - 1) (1-2) T1() e((2m-1)2m/4N]t- (16) m= (2m -1)2m For t > 0, the series is easily seen to be convergent, and, as t -a c, %:goes to zero. For t = 0, we must show that this series converges to 1. Let o0, n be a partial sum of the first n terms; then, by the recurrence relation (4m - 1) (1 - r2)T2'm_2(r)/ (2m - 1)2m = P2m_2(r) - P2m(r), we have Ro, n = 1 - P2n(r).- By using an in- tegral expression of Pn, i.e., Pn(z) = (1/r) fo I{Z + 2 - 1 COS t} n dt, we can show that, for IrI < 1 P2,(r) -O 0 as n co. For IP2n(r) I < - jr + Vr2 - 1 COS t 1 2? dt = -1fT {r2 + (1 -r2) Cos2 t dt 0 (n co). Furthermore, from (16), co = = _ {P2J(r) - P2j+2 (r)}eI[(2+l)(2i+2)/4 t. (17) Consequently, it can be seen that, from equations (8) and (17), f(l, t) + Qt + f(O, t) = 1, as it should be. The processes of the change in the distribution of the unfixed classes when the population starts from p = 0.5 and p = 0.1 are illustrated in Figures 1 and 2, respec- tively. In Figure 1 it will be seen that after 2N generations the distribution curve becomes almost flat, and the genes are still unfixed in about 50 per cent of the cases. In Figure 2, the initial gene frequency is assumed to be 10 per cent, and it takes 4N or 5N generations before the distribution curve becomes practically flat. By that time, however, the genes are fixed in more than 90 per cent of the cases, and the simplest asymptotic formula Ce- (1/2N)t may not be useful as in the case of p = 0.5. 148 -P-Roc. N. A. S. GENETICS: MOTOO KIMURA FIGS. 1-2.-The processes of the change in the probability distribution of heterallelic classes, due to random sampling of gametes in reproduction. It is assumed that the population starts from the gene frequency 0.5 in Fig. 1 (left) and 0.1 in Fig. 2 (right). t = time in genera- tion; N = effective size of the population; abscissa is gene frequency; ordinate is probability density. The probability of heterozygosis is calculated by equation (15): fo12x(1-x~~x~t~dx =(2i+ l)TH. = O 2X(1- X)+O(X, t)di= E (i+1) i i (1-2p) X (1 -Z2)T,._.(z) e-'i(i + 1)/4N]t dZ. By virtue of equation (14) (put m = 0), the last integral is 0 except for i = 1. Hence Hg= pq 1 4 (2)t 2pqe-(l/2lt = Hoe-(l/2N)t, (18)2 3 showing that the heterozygosis decreases exactly at the rate of 1/(2N) per generation. This is readily confirmed by a simple calculation: Let p be the frequency of A in the population, where the frequency of the heterozygotes is 2p(l - p). The amount of heterozygosis to be expected after one generation of random sampling of the gametes is E{2 ( + 5p) (1 -P- 6P)} 2p(1 -p) -2E(ap)2= 2p(l - p) - 2 =(-22p(1-p), as was to be shown. 149VOL. 41) 1955 GENETICS: SAX AND KING The author expresses his appreciation to Dr. James F. Crow for valuable help during the course of this work. Thanks are also due to Dr. E. R. Immel for his helpful suggestions. * Contribution No. 84 of the National Institute of Genetics, Mishima-shi, Japan. t Contribution No. 570. This work was supported by a grant from the University Research Committee from funds supplied by the Wisconsin Alumni Research Foundation. 1 R. A. Fisher, Proc. Roy. Soc. Edinburgh, 42,321-341, 1922. 2 S. Wright, Genetics, 16, 97-159, 1931. 3 R. A. Fisher, The Genetical Theory of Natural Selection (Oxford, 1930). 4M. Kimura, Genetics, 39, 280-295, 1954. 6 S. Goldberg, Ph.D. thesis, Cornell University, 1950. 6 A. Robertson, Genetics, 37, 189-207, 1952. 7M. Morse and H. Feshbach, Methods of Theoretical Physics (New York, 1953). 8 S. Wright, these PROCEEDINGS, 31,382-389, 1945. 9 W. Feller, Proc. Second Berkeley Symposium on Math. Stat. and Prob., Univ. of California, pp. 227-246, 1951. AN X-RAY ANALYSIS OF CHROMOSOME DUPLICATION* BY KARL SAX AND EDWARD D. KING ARNOLD ARBORETUM AND BUSSEY INSTITUTION, HARVARD UNIVERSITY, JAMAICA PLAIN, MASSACHUSETTS Communicated January 5, 1955 The induction of half-chromatid aberrations by X-rays provides evidence re- garding the multiple nature of the chromosome and may be of considerable signifi- cance in the analysis of chromosome duplication. The occurrence of half-chro- matid exchanges was first described by Swanson' in 1943. These aberrations, found in pollen-tube mitoses of Tradescantia, were of sporadic occurrence. Re- cently Crouse2 has found that practically all the X-ray-induced aberrations induced at the first meiotic metaphase stage in Lilium involve half-chromatid breaks and exchanges. Breaks in one of the half-chromatids of each of two chromatids are followed by reciprocal translocation to produce a chromatid bridge at anaphase. The close association of the coiled half-chromatids prevents the terminal separation of the two anaphase chromosomes (cf. Fig. 1). Occasional half-chromatid aberrations have been found at anaphase in the divi- sion of the microspore nucleus of Tradescantia following X-irradiation. The dosages commonly used (50-150 r) produced considerable stickiness of the chromo- somes, so that accurate analysis of chromosomal aberrations could not be made until about 6 hours after raying. By reducing the dosage to 25 r and keeping the inflorescences at 30 C. during irradiation, it was possible to obtain clear figures of mitoses and induced aberrations as early as 3 hours after raying. Irradiation at 30 C. was done to compensate for the low dosage, since it has been shown that the chromosome aberration frequency can be infeased by raying at low temperatures.3 Chromosomal aberrations induced at 4-6 hours before anaphasA consisted almost entirely of half-chromatid bridges at anaphase. Presumably, these chromosomes were irradiated at prometaphase or metaphase. At these stages the two sister 150 PROC. N. A. S.