5.5: Випадковий генетичний дрейф

До цих пір наші спрощені генетичні моделі передбачали нескінченну популяцію, нехтуючи стохастичними ефектами. Тут розглядається кінцева популяція: результуючі стохастичні ефекти на алельні частоти називають випадковим генетичним дрейфом. Найпростіша генетична модель, що включає випадковий генетичний дрейф, передбачає популяції\(N\) особин фіксованого розміру та моделює еволюцію діалельного гаплоїдного генетичного локусу.

Існує дві широко використовувані генетичні моделі для скінченних популяцій: модель РайтФішера та модель Морана. Модель Райта-Фішера найбільш схожа на нашу модель дискретного покоління з нескінченною популяцією. У моделі Райта-Фішера\(N\) дорослі особини виділяють дуже велику кількість гамет в генофонд, а наступне покоління формується з\(N\) випадкових гамет, самостійно обраних з генофонду. Модель Морана використовує інший підхід. У моделі Морана один крок еволюції складається з однієї випадкової індивіди в популяції, що відтворюється, а інша випадкова особина вмирає, при цьому популяція завжди підтримується в постійному розмірі\(N\). Оскільки дві випадкові події відбуваються в кожному поколінні в моделі Морана, а\(N\) випадкові події відбуваються в кожному поколінні в моделі Райта-Фішера, загальна кількість етапів\(N / 2\) еволюції в моделі Морана порівнянна, але не зовсім ідентична, з одним дискретним поколінням в Модель Райта-Фішера. Однак було показано, що дві моделі стають ідентичними в межі великих\(N\). Для наших цілей модель Морана математично більш зрозуміла, і ми приймаємо її тут.

Ми розробляємо нашу модель аналогічно стохастичній моделі зростання популяції, отриманої в\(\S 3.1\). \(n\)Позначимо кількість\(A\) -алелей в популяції,\(N-n\) а кількість\(a\) -алелей. При\(n=0,1,2, \ldots, N\) дискретній випадковій величині функція маси ймовірності\(p_{n}(g)\) позначає ймовірність наявності\(n A\) -алелей на етапі еволюції\(g .\) З\(n\) А-алелями в популяції розміром \(N\), ймовірність особи, яка відтворює або вмирає, задається\(s_{n}=n / N\); відповідна ймовірність для\(a\) є\(1-s_{n}\).\(A\) Існує три способи отримання популяції\(n A\) -алелей на етапі еволюції\(g+1\). По-перше, були\(n A\) -алелі на етапі еволюції\(g\), і індивіда, що відтворює, несла той самий алель, що і індивід, що вмирає. По-друге, на етапі еволюції\(n-1\) були -алелі\(g\), і індивід, що відтворює,\(A\) і індивід, що вмирає\(a\). І по-третє, були\(n+1 A\) -алелі на етапі еволюції\(g\), і індивід, що відтворює має,\(a\) і індивід вмирає\(A\). Множення ймовірностей та підсумовування трьох випадків призводить до

\[\begin{align} \nonumber p_{n}(g+1)=\left(s_{n}^{2}+\left(1-s_{n}\right)^{2}\right) p_{n}(g)+s_{n-1}\left(1-s_{n-1}\right) p_{n-1}(g) & \\[4pt] &+s_{n+1}\left(1-s_{n+1}\right) p_{n+1}(g) \end{align} \nonumber \]

Зауважте, що це рівняння дійсне для\(0<n<N\), і що рівняння на границях, що представляють ймовірності того, що один з алелей фіксований - є

\[\begin{align} p_{0}(g+1) &=p_{0}(g)+s_{1}\left(1-s_{1}\right) p_{1}(g) \\[4pt] p_{N}(g+1) &=p_{N}(g)+s_{N-1}\left(1-s_{N-1}\right) p_{N-1}(g) \end{align} \nonumber \]

Межі називаються абсорбуючими і ймовірність фіксації алеля монотонно зростає з кожним народженням і смертю. Як тільки ймовірність фіксації алеля дорівнює одиниці, подальших змін частот алелів не відбувається.

Проілюструємо рішення\((5.5.1)\) - (5.5.3) на рис. \(5.3\)для малої чисельності населення\(N=20\), і де число\(A\) -алелей точно відомо в поколінні заснування, з або (а)\(p_{10}(0)=1\), або; (б)\(p_{13}(0)=1\). Побудовано масову густину ймовірності кількості\(A\) -алелей на кожному етапі\(N\) еволюції, аж до\(7 N\) кроків, що відповідають приблизно чотирнадцяти дискретним поколінням еволюції в моделі Райта-Фішера. Зверніть увагу, як розподіл ймовірностей розсіюється від його початкового значення, і як ймовірності в кінцевому підсумку концентруються на кордоні, причому обидва\(p_{0}\) і\(p_{20}\) монотонно збільшуються. Насправді після великої кількості поколінь\(p_{0}\) наближається до початкової частоти\(a\) -аллеля і\(p_{20}\) наближається до початкової частоти\(A\) -аллеля (на малюнках не показана).

Малюнок 5.3:\(p_{n}\) проти\(n\) с\(N=20\). Еволюції етапи побудовані відповідають\(g=\)\(0, N, 2 N, \ldots, 7 N\). (а) Початкова кількість\(A\) осіб є\(n=10\). (б) Початкова кількість\(A\) осіб\(n=13\) є. і зауважте, що\(S(x)=x .\) Аналогічно, ми маємо\(s_{n+1}=S(x+\Delta x)\) і\(s_{n-1}=S(x-\Delta x)\), де\(\Delta x=1 / N .\) Тоді, з\(\Delta t=2 / N,(5.5.4)\) перетворюється в<\ p>

\[\begin{align} \nonumber P(x, t&+\Delta t)-P(x, t)=S(x+\Delta x)(1-S(x+\Delta x)) P(x+\Delta x, t) \\[4pt] &-2 S(x)(1-S(x)) P(x, t)+S(x-\Delta x)(1-S(x-\Delta x)) P(x-\Delta x, t) \end{align} \nonumber \]

Для подальшого спрощення використано відоме центрально-різницеве наближення до другої похідної функції\(f(x)\),

\[f^{\prime \prime}(x)=\frac{f(x+\Delta x)-2 f(x)+f(x-\Delta x)}{\Delta x^{2}}+\mathrm{O}\left(\Delta x^{2}\right) \nonumber \]

і розпізнати праву частину (5.5.7) чисельником другої похідної. З\(\Delta x=\Delta t / 2=1 / N \rightarrow 0\), і\(S(x)=x\), виводимо до порядку відведення в рівнянні\(1 / N\) з частинними похідними

\[\frac{\partial P(x, t)}{\partial t}=\frac{1}{2} \frac{\partial^{2}}{\partial x^{2}}(V(x) P(x, t)) \nonumber \]

із

\[V(x)=\frac{x(1-x)}{N} \nonumber \]

Для інтерпретації значення функції\(V(x)\) ми будемо використовувати наступний результат теорії ймовірностей. Для\(n\) самостійних випробувань кожен з ймовірністю успіху\(p\) і ймовірністю невдачі\(1-p\) число успіхів, що позначається\(X\), є біноміальною випадковою величиною з параметрами\((n, p)\). Загальновідомими результатами\(E[\ldots]\) є\(E[X]=n p\) і\(\operatorname{Var}[X]=n p(1-p)\), де очікуване значення, і\(\operatorname{Var}[\ldots]\) є дисперсія.

Тепер кількість\(A\) -алелей, обраних при формуванні наступного покоління Райта-Фішера, є біноміальною випадковою величиною\(n^{\prime}\) з параметрами\((N, n / N)\). Тому\(E\left[n^{\prime}\right]=n\) і\(\operatorname{Var}\left[n^{\prime}\right]=n(1-n / N)\). З\(x^{\prime}=n^{\prime} / N\) і\(x=n / N\), у нас\(E\left[x^{\prime}\right]=x\), і\(\operatorname{Var}\left[x^{\prime}\right]=x(1-x) / N\). Тому функцію\(V(x)\) можна інтерпретувати як дисперсію\(x\) над одним поколінням Райта-Фішера.

Хоча (5.5.9) і (5.5.10) явно залежать від чисельності населення\(N\), чисельність населення може бути усунена простою зміною змінних. Якщо ми дозволимо

\[\tau=t / N \text {, } \nonumber \]

тоді диференціальне рівняння (5.5.9) перетворюється на

\[\frac{\partial P(x, \tau)}{\partial \tau}=\frac{1}{2} \frac{\partial^{2}}{\partial x^{2}}(x(1-x) P(x, \tau)) \nonumber \]

незалежний від\(N\). Зміна змінних, задана (5.5.11), просто стверджує, що подвоєння чисельності популяції подовжить часову шкалу еволюції на відповідний коефіцієнт у два рази. Пам'ятайте, що тут ми вже працюємо під припущенням, що чисельність населення велика.

Рівняння (5.5.12) являє собою дифузійно-подібне рівняння для функції розподілу ймовірностей\(P(x, \tau)\). Дифузійне наближення для вивчення генетичного дрейфу було вперше введено в популяційну генетику її засновниками Фішером і Райтом, а пізніше широко розроблено в\(1950^{\prime}\) s японським біологом Motoo Kimura. Тут наш аналіз цього рівняння спирається на більш свіжу роботу Маккейна і Ваксмана (2007), які показали, як побудувати аналітичне рішення (5.5.12) на межах\(x\). Бо\(0<x<1\), це наводить на думку з числових рішень, показаних на рис. \(5.3\)що розподіл ймовірностей може асимптотично стати незалежним від\(x\). Відповідно, шукаємо асимптотичне рішення для\((5.5.12)\) внутрішніх значень\(x\) задовольняють\(P(x, \tau)=P(\tau)\). Рівняння (5.5.12) стає

\[\begin{aligned} \frac{d P(\tau)}{d \tau} &=\frac{1}{2}(x(1-x))^{\prime \prime} P(\tau) \\[4pt] &=-P(\tau) \end{aligned} \nonumber \]

з асимптотичним розчином

\[P(\tau)=c e^{-\tau} \nonumber \]

і ми спостерігаємо, що сумарна ймовірність над внутрішніми значеннями\(x\) розпадається в геометричній прогресії.

Щоб зрозуміти, як поводиться рішення на кордоні\(x\), потрібні граничні умови для\(P(x, \tau)\). Насправді, оскільки\(P(x, \tau)\) є одниною на кордоні\(x\), відповідні граничні умови не можуть бути отримані безпосередньо з різницевих рівнянь, заданих (5.5.2 і 5.5.3). Швидше за все граничні умови найлегше отримати при першому переказуванні (5.5.12) у вигляді рівняння неперервності.

Давайте\(j(x, \tau)\) позначимо так званий струм ймовірності. У невеликій області розміру,\(\Delta x\) що лежить в інтервалі\((x, x+\Delta x)\), часова швидкість зміни ймовірності\(P(x, \tau) \Delta x\) обумовлена потоком ймовірності в цю область і з неї. При відповідному визначенні ймовірності струму ми маємо в загальному

\[\frac{\partial}{\partial \tau}(P(x, \tau) \Delta x)=j(x, t)-j(x+\Delta x) \nonumber \]

або як\(\Delta x \rightarrow 0\),

\[\frac{\partial P(x, \tau)}{\partial \tau}+\frac{\partial j(x, \tau)}{\partial x}=0 \nonumber \]

яка є звичайною формою рівняння неперервності. Ідентифікація (5.5.15) з (5.5.12) показує, що ймовірність струму нашої задачі задається

\[j(x, \tau)=-\frac{1}{2} \frac{\partial}{\partial x}(x(1-x) P(x, \tau)) . \nonumber \]

Тепер, оскільки повна ймовірність є єдністю; тобто

\[\int_{0}^{1} P(x, \tau) d x=1 \nonumber \]

ймовірність не може протікати в межах або за межі\(x\), і тому ми повинні мати

\[j(0, \tau)=0, \quad j(1, \tau)=0, \nonumber \]

які є відповідними граничними умовами для (5.5.12).

Ми можемо шукати стаціонарні рішення (5.5.12). Використання рівняння неперервності (5.5.15) разом з граничними умовами (5.5.18) показує, що стаціонарне рішення має нульовий струм ймовірності; тобто\(j(x)=0\). Інтеграція\(j(x)\) використання (5.5.16) призводить до

\[x(1-x) P(x)=c_{1} \nonumber \]

де\(c_{1}\) - константа інтеграції. Легко очевидним рішенням є\(P(x)=\)\(c_{1} /[x(1-x)]\), але є також два менш очевидних рішення. Функціям розподілу ймовірностей дозволено містити сингулярні розв'язки, що відповідають дельта-функціям Дірака, і тут розглядається можливість існування дельта-функцій Дірака на границях. Помітивши, що обидва\(x \delta(x)=0\) і\((1-x) \delta(1-x)=0\), ми можемо побачити, що загальне рішення для\(P(x)\) може бути записано як

\[P(x)=\frac{c_{1}}{x(1-x)}+c_{2} \delta(x)+c_{3} \delta(1-x) . \nonumber \]

\(P(x)\)Вимога задовольнити (5.5.17) призводить до\(c_{1}=0\) і\(c_{2}+c_{3}=1\). Щоб визначити залишилася вільну константу, ми можемо обчислити середню частоту\(A\) -аллеля в популяції. З рівняння\((5.5.15)\) неперервності ми маємо

\[\begin{aligned} \frac{\partial}{\partial \tau} \int_{0}^{1} x P(x, \tau) d x &=-\int_{0}^{1} x \frac{\partial j(x, \tau)}{\partial x} d x \\[4pt] &=\int_{0}^{1} j(x, \tau) d x \\[4pt] &=0 \end{aligned} \nonumber \]

де перший інтеграл з правого боку був зроблений частинами з використанням граничних умов, заданих (5.5.18), а другий інтеграл був зроблений за допомогою (5.5.16) і зникнення\(x(1-x) P(x)\) на кордоні. Отже, середня частота\(A\) алеля є постійною, як можна було б очікувати для ненаправленого випадкового генетичного дрейфу - і ми можемо припустити, що його початкове значення є\(p\). Тому ми отримуємо для нашого стаціонарного рішення

\[P(x)=(1-p) \delta(x)+p \delta(1-x) \nonumber \]

Отже, можлива ймовірність фіксації\(A\) -алеля просто дорівнює його початковій частоті. Наприклад, припустимо, що в межах популяції однорідних\(N\) для однорідної для\(a\) однієї нейтральної мутації відбувається так, що одна особина тепер несе\(A\) -аллель. Яка ймовірність того, що\(A\) -аллель врешті-решт стане фіксованим? Наш результат дасть ймовірність\(1 / N\), яка є початковою частотою\(A\) -аллеля. Інтуїтивно, після того, як пройшла достатня кількість поколінь, всі живі особини повинні бути нащадками від єдиного родового індивіда, що живе в той момент, коли сталася єдина мутація. Імовірність того, що ця одинична особина несла\(A\) -аллель, справедлива\(1 / N\).

Відзначимо далі, що Кімура першим знайшов аналітичне рішення рівняння дифузії (5.5.12). Метод розв'язання з використанням перетворень Фур'є можна знайти в додатку Маккейна і Ваксмана (2007). Окрім використання дельта-функцій Дірака, цим авторам потрібні крокові функції Хевісайда, функції Бесселя, сферичні функції Бесселя, гіпергеометричні функції, поліноми Лежандра та поліноми Гегенбауера. Отриманий розчин має форму

\[P(x, \tau)=\Pi_{0}(\tau) \delta(x)+\Pi_{1}(\tau) \delta(1-x)+f(x, \tau) \nonumber \]

Якщо кількість\(A\) -алелей відомо в початковий момент, а частота є\(p\), то

\[P(x, 0)=\delta(x-p) \nonumber \]

і\(\Pi_{0}(0)=\Pi_{1}(0)=0 ; f(x, 0)=\delta(x-p)\). Як\(\tau \rightarrow \infty\), у нас є\(f(x, \tau) \rightarrow 0\);\(\Pi_{0}(\tau) \rightarrow 1-p\) і\(\Pi_{1}(\tau) \rightarrow p\)

Ми можемо легко продемонструвати хоча б одне точне рішення форми (5.5.22). Для початкового рівномірного розподілу ймовірностей, заданого

\[P(x, 0)=1, \nonumber \]

розподіл ймовірностей у внутрішніх точках залишається рівномірним і задається\((5.5.13)\) з\(c=1\). Якщо взяти форму розв'язку, задану (5.5.22) разом з вимогою одиничної ймовірності, заданої (5.5.17), ми можемо отримати точний результат

\[P(x, \tau)=\frac{1}{2}\left(1-e^{-\tau}\right)(\delta(x)+\delta(1-x))+e^{-\tau} \nonumber \]