3.2: Asintóticos de grandes poblaciones iniciales

Nuestro objetivo aquí es resolver una expansión de la distribución en poderes de\(1 / N_{0}\) a órdenes principales; observe que\(1 / N_{0}\) es pequeño si\(N_{0}\) es grande. Al orden cero, es decir en el límite\(N_{0} \rightarrow \infty\), mostraremos que se recupera el modelo determinista de crecimiento poblacional. A primer orden en\(1 / N_{0}\), mostraremos que la distribución de probabilidad es normal. Este último resultado se verá como consecuencia del conocido Teorema del Límite Central en la teoría de la probabilidad.

Desarrollamos nuestra expansión trabajando directamente con la ecuación diferencial para\(p_{N}(t)\). Ahora, cuando el tamaño de la población\(N\) es una variable aleatoria discreta (tomando solo los valores enteros no negativos de\(0,1,2, \ldots), p_{N}(t)\) es la función de masa de probabilidad para\(N .\) Si\(N_{0}\) es grande, entonces la naturaleza discreta de\(N\) es intrascendente, y es preferible trabajar con una variable aleatoria continua y su función de densidad de probabilidad. En consecuencia, definimos la variable aleatoria\(x=N / N_{0}\), y la tratamos\(x\) como una variable aleatoria continua, con\(0 \leq x<\infty\). Ahora,\(p_{N}(t)\) es la probabilidad de que la población sea de tamaño\(N\) en el momento\(t\), y la función de densidad de probabilidad de\(x, P(x, t)\), se define tal que\(\int_{a}^{b} P(x, t) d x\) es la probabilidad de que\(a \leq x \leq b\). La relación entre\(p\) y se\(P\) puede determinar considerando cómo aproximar una distribución discreta de probabilidad por una distribución continua, es decir, definiendo\(P\) tal que

\[\begin{aligned} p_{N}(t) &=\int_{\left(N-\frac{1}{2}\right) / N_{0}}^{\left(N+\frac{1}{2}\right) / N_{0}} P(x, t) d x \\[4pt] &=P\left(N / N_{0}, t\right) / N_{0} \end{aligned} \nonumber \]

donde la última igualdad se vuelve exacta como\(N_{0} \rightarrow \infty\). Por lo tanto, la definición apropiada para\(P(x, t)\) viene dada por

\[P(x, t)=N_{0} p_{N}(t), \quad x=N / N_{0} \nonumber \]

que satisface

\[\begin{aligned} \int_{0}^{\infty} P(x, t) d x &=\sum_{N=0}^{\infty} P\left(N / N_{0}, t\right)\left(1 / N_{0}\right) \\[4pt] &=\sum_{N=0}^{\infty} p_{N}(t) \\[4pt] &=1 \end{aligned} \nonumber \]

la primera igualdad (exacta sólo cuando\(N_{0} \rightarrow \infty\)) siendo una suma de Reimann aproximación de la integral.

Ahora transformamos el conjunto infinito de ecuaciones diferenciales ordinarias (3.1.3) para\(p_{N}(t)\) en una única ecuación diferencial parcial para\(P(x, t)\). Multiplicamos (3.1.3) por\(N_{0}\) y sustituimos\(N=N_{0} x, p_{N}(t)=P(x, t) / N_{0}\), y\(p_{N-1}(t)=P\left(x-\frac{1}{N_{0}}, t\right) / N_{0}\) para obtener

\[\frac{\partial P(x, t)}{\partial t}=b\left[\left(N_{0} x-1\right) P\left(x-\frac{1}{N_{0}}, t\right)-N_{0} x P(x, t)\right] \nonumber \]

A continuación, las series de Taylor se expanden\(P\left(x-1 / N_{0}, t\right)\) alrededor\(x\), tratando\(1 / N_{0}\) como un pequeño parámetro. Es decir, hacemos uso de

\[\begin{aligned} P\left(x-\frac{1}{N_{0}}, t\right) &=P(x, t)-\frac{1}{N_{0}} P_{x}(x, t)+\frac{1}{2 N_{0}^{2}} P_{x x}(x, t)-\ldots \\[4pt] &=\sum_{n=0}^{\infty} \frac{(-1)^{n}}{n ! N_{0}^{n}} \frac{\partial^{n} P}{\partial x^{n}} . \end{aligned} \nonumber \]

Los dos términos principales proporcionales a\(N_{0}\) en el lado derecho de (3.2.2) cancelan exactamente, y si agrupamos los términos restantes en poderes de\(1 / N_{0}\), obtenemos para los tres primeros términos principales en la expansión

\[\begin{align} \nonumber P_{t} &=-b\left[\left(x P_{x}+P\right)-\frac{1}{N_{0}}\left(\frac{x P_{x x}}{2 !}+\frac{P_{x}}{1 !}\right)+\frac{1}{N_{0}^{2}}\left(\frac{x P_{x x x}}{3 !}+\frac{P_{x x}}{2 !}\right)-\ldots\right] \\[4pt] &=-b\left[(x P)_{x}-\frac{1}{N_{0} 2 !}(x P)_{x x}+\frac{1}{N_{0}^{2} 3 !}(x P)_{x x x}-\ldots\right] \end{align} \nonumber \]

y los términos de orden superior se pueden obtener siguiendo el patrón evidente.

La ecuación (3.2.3) puede analizarse adicionalmente mediante una expansión de perturbación de la función de densidad de probabilidad en potencias de\(1 / N_{0}\):

\[P(x, t)=P^{(0)}(x, t)+\frac{1}{N_{0}} P^{(1)}(x, t)+\frac{1}{N_{0}^{2}} P^{(2)}(x, t)+\ldots \nonumber \]

Aquí, las funciones desconocidas\(P^{(0)}(x, t), P^{(1)}(x, t), P^{(2)}(x, t)\), etc. se van a determinar sustituyendo la expansión (3.2.4) en\((3.2.3)\) e igualando coeficientes de potencias de\(1 / N_{0}\). Así obtenemos para los coeficientes de\(\left(1 / N_{0}\right)^{0}\) y\(\left(1 / N_{0}\right)^{1}\),

\[\begin{align} P_{t}^{(0)} &=-b\left(x P^{(0)}\right)_{x}^{\prime} \\[4pt] P_{t}^{(1)} &=-b\left[\left(x P^{(1)}\right)_{x}-\frac{1}{2}\left(x P^{(0)}\right)_{x x}\right] \end{align} \nonumber \]

Derivación del modelo determinista

El término de orden cero en la expansión de perturbación (3.2.4),

\[P^{(0)}(x, t)=\lim _{N_{0} \rightarrow \infty} P(x, t) \nonumber \]

satisface (3.2.5). La ecuación (3.2.5) es una ecuación diferencial parcial lineal de primer orden con un coeficiente variable. Una forma de resolver esta ecuación es probar el ansatz

\[P^{(0)}(x, t)=h(t) f(r), \quad r=r(x, t), \nonumber \]

junto con la condición inicial

\[P^{(0)}(x, 0)=f(x) \nonumber \]

ya que at\(t=0\), se supone que la distribución de probabilidad es una función conocida. Esta condición inicial impone otras restricciones útiles a las funciones\(h(t)\) y\(r(x, t)\):

\[h(0)=1, \quad r(x, 0)=x \nonumber \]

Las derivadas parciales de (3.2.8) son

\[P_{t}^{(0)}=h^{\prime} f+h r_{t} f^{\prime}, \quad P_{x}^{(0)}=h r_{x} f^{\prime}, \nonumber \]

que después de la sustitución en (3.2.5) da como resultado

\[\left(h^{\prime}+b h\right) f+\left(r_{t}+b x r_{x}\right) h f^{\prime}=0 \nonumber \]

Esta ecuación puede ser satisfecha para cualquier\(f\) siempre que\(h(t)\) y\(r(x, t)\) satisfaga

\[h^{\prime}+b h=0, \quad r_{t}+b x r_{x}=0 \nonumber \]

La primera ecuación para\(h(t)\), junto con la condición inicial\(h(0)=1\), se resuelve fácilmente para producir

\[h(t)=e^{-b t} . \nonumber \]

Para determinar una solución para\(r(x, t)\), intentamos la técnica de separación de variables. Escribimos\(r(x, t)=X(x) T(t)\), y tras la sustitución en la ecuación diferencial para\(r(x, t)\), obtenemos

\[X T^{\prime}+b x X^{\prime} T=0 \nonumber \]

y división por XT y rendimientos de separación

\[\frac{T^{\prime}}{T}=-b x \frac{X^{\prime}}{X} \nonumber \]

Dado que el lado izquierdo es independiente de\(x\), y el lado derecho es independiente de\(t\), tanto el lado izquierdo como el lado derecho deben ser constantes, independientes de ambos\(x\) y\(t .\) Ahora, nuestra condición inicial es \(r(x, 0)=x\), para que\(X(x) T(0)=x .\) Sin pérdida de generalidad, podamos tomar\(T(0)=1\), para que\(X(x)=x .\) El lado derecho de (3.2.16) sea por lo tanto igual a la constante\(-b\), y obtengamos la ecuación diferencial y la condición inicial

\[T^{\prime}+b T=0, \quad T(0)=1 \nonumber \]

que podemos resolver para ceder\(T(t)=e^{-b t}\). Por lo tanto, nuestra solución para\(r(x, t)\) es

\[r(x, t)=x e^{-b t} \nonumber \]

Al poner nuestras soluciones (3.2.14) y (3.2.18) juntas en nuestro ansatz (3.2.8), hemos obtenido la solución general a la pde:

\[P^{(0)}(x, t)=e^{-b t} f\left(x e^{-b t}\right) \nonumber \]

Para determinar\(f\), aplicamos la condición inicial de la función de masa de probabilidad,\(p_{N}(0)=\delta_{N, N_{0}}\). De (3.10), la condición inicial correspondiente en la función de distribución de probabilidad es

\[P(x, 0)= \begin{cases}N_{0} & \text { if } 1-\frac{1}{2 N_{0}} \leq x \leq 1+\frac{1}{2 N_{0}} \\[4pt] 0 & \text { otherwise }\end{cases} \nonumber \]

En el límite\(N_{0} \rightarrow \infty, P(x, 0) \rightarrow P^{(0)}(x, 0)=\delta(x-1)\), donde\(\delta(x-1)\) está la delta-función Dirac, centrada alrededor de 1. La función delta es ampliamente utilizada en la física cuántica y fue introducida por Dirac para ese propósito. Ahora encuentra muchos usos en las matemáticas aplicadas. Se puede definir exigiendo que, para cualquier función\(g(x)\),

\[\int_{-\infty}^{+\infty} g(x) \delta(x) d x=g(0) \nonumber \]

La visión habitual de la función delta\(\delta(x-a)\) es que es cero en todas partes excepto\(x=a\) en la que es infinita, y su integral es una. No es realmente una función, sino que es lo que los matemáticos llaman una distribución.

Ahora, ya que\(P^{(0)}(x, 0)=f(x)=\delta(x-1)\), nuestra solución se convierte en

\[P^{(0)}(x, t)=e^{-b t} \delta\left(x e^{-b t}-1\right) \nonumber \]

Esto puede ser reescrito señalando que (dejar\(y=a x-c\)),

\[\begin{aligned} \int_{-\infty}^{+\infty} g(x) \delta(a x-c) d x &=\frac{1}{a} \int_{-\infty}^{+\infty} g((y+c) / a) \delta(y) d y \\[4pt] &=\frac{1}{a} g(c / a) \end{aligned} \nonumber \]

cediendo la identidad

\[\delta(a x-c)=\frac{1}{a} \delta\left(x-\frac{c}{a}\right)\nonumber \]

A partir de esto, podemos reescribir nuestra solución\((3.2.22)\) de la forma más intuitiva

\[P^{(0)}(x, t)=\delta\left(x-e^{b t}\right) \nonumber \]

Usando\((3.2.23)\), el valor esperado de orden cero de\(x\) es

\[\begin{aligned} \left\langle x_{0}\right\rangle &=\int_{0}^{\infty} x P^{(0)}(x, t) d x \\[4pt] &=\int_{0}^{\infty} x \delta\left(x-e^{b t}\right) d x \\[4pt] &=e^{b t} \end{aligned} \nonumber \]

mientras que la varianza de orden cero es

\[\begin{aligned} \sigma_{x_{0}}^{2} &=\left\langle x_{0}^{2}\right\rangle-\left\langle x_{0}\right\rangle^{2} \\[4pt] &=\int_{0}^{\infty} x^{2} P^{(0)}(x, t) d x-e^{2 b t} \\[4pt] &=\int_{0}^{\infty} x^{2} \delta\left(x-e^{b t}\right) d x-e^{2 b t} \\[4pt] &=e^{2 b t}-e^{2 b t} \\[4pt] &=0 \end{aligned} \nonumber \]

Así, en el límite infinito de población, la variable aleatoria\(x\) tiene varianza cero, y por lo tanto ya no es aleatoria, sino que sigue\(x=e^{b t}\) determinísticamente. Decimos que la distribución de probabilidad de\(x\) se vuelve aguda en el límite de grandes tamaños de población. El principio general de modelar poblaciones grandes determinísticamente puede simplificar los modelos matemáticos cuando los efectos estocásticos no son importantes.

Derivación de la distribución de probabilidad normal

Consideramos ahora el término de primer orden en la expansión de perturbación (3.2.4), que satisface (3.2.6). No sabemos resolver (3.2.6) directamente, por lo que intentaremos encontrar una solución siguiendo una ruta más sinuosa. Primero, se procede calculando los momentos de la distribución de probabilidad. Tenemos

\[\begin{aligned} \left\langle x^{n}\right\rangle &=\int_{0}^{\infty} x^{n} P(x, t) d x \\[4pt] &=\int_{0}^{\infty} x^{n} P^{(0)}(x, t) d x+\frac{1}{N_{0}} \int_{0}^{\infty} x^{n} P^{(1)}(x, t) d x+\ldots \\[4pt] &=\left\langle x_{0}^{n}\right\rangle+\frac{1}{N_{0}}\left\langle x_{1}^{n}\right\rangle+\ldots, \end{aligned} \nonumber \]

donde define la última igualdad\(\left\langle x_{0}^{n}\right\rangle\), etc. Ahora, usando (3.2.23),

\[\begin{align} \nonumber \left\langle x_{0}^{n}\right\rangle &=\int_{0}^{\infty} x^{n} P^{(0)}(x, t) d x \\[4pt] \nonumber &=\int_{0}^{\infty} x^{n} \delta\left(x-e^{b t}\right) d x \\[4pt] &=e^{n b t} \end{align} \nonumber \]

Para determinar\(\left\langle x_{1}^{n}\right\rangle\), utilizamos\((3.2.6)\). Multiplicando\(x^{n}\) e integrando, tenemos

\[\frac{d\left\langle x_{1}^{n}\right\rangle}{d t}=-b\left[\int_{0}^{\infty} x^{n}\left(x P^{(1)}\right)_{x} d x-\frac{1}{2} \int_{0}^{\infty} x^{n}\left(x P^{(0)}\right)_{x x} d x\right] \nonumber \]

Nos integramos por partes para eliminar las derivadas de\(x P\), asumiendo que\(x P\) y todas sus derivadas desaparecen en los límites de la integración, donde\(x\) es igual a cero o infinito. Tenemos

\[\begin{aligned} \int_{0}^{\infty} x^{n}\left(x P^{(1)}\right)_{x} d x &=-\int_{0}^{\infty} n x^{n} P^{(1)} d x \\[4pt] &=-n\left\langle x_{1}^{n}\right\rangle \end{aligned} \nonumber \]

y

\[\begin{aligned} \frac{1}{2} \int_{0}^{\infty} x^{n}\left(x P^{(0)}\right)_{x x} d x &=-\frac{n}{2} \int_{0}^{\infty} x^{n-1}\left(x P^{(0)}\right)_{x} d x \\[4pt] &=\frac{n(n-1)}{2} \int_{0}^{\infty} x^{n-1} P^{(0)} d x \\[4pt] &=\frac{n(n-1)}{2}\left\langle x_{0}^{n-1}\right\rangle \end{aligned} \nonumber \]

Por lo tanto, después de la integración por partes, (3.2.25) se convierte en

\[\frac{d\left\langle x_{1}^{n}\right\rangle}{d t}=-b\left[n\left\langle x_{1}^{n}\right\rangle+\frac{n(n-1)}{2}\left\langle x_{0}^{n-1}\right\rangle\right] \nonumber \]

La ecuación (3.2.26) es una ecuación diferencial lineal no homogénea de primer orden y se puede resolver usando un factor integrador (ver (3.1.10) y la discusión anterior). Resolviendo esta ecuación diferencial usando (3.2.24) y la condición inicial\(\left\langle x_{1}^{n}\right\rangle(0)=0\), obtenemos

\[\left\langle x_{1}^{n}\right\rangle=\frac{n(n-1)}{2} e^{n b t}\left(1-e^{-b t}\right) \nonumber \]

La función de distribución de probabilidad, exacta al orden\(1 / N_{0}\), puede obtenerse haciendo uso de la llamada función de generación de momentos\(\Psi(s)\), definida como

\[\begin{aligned} \Psi(s) &=\left\langle e^{s x}\right\rangle \\[4pt] &=1+s\langle x\rangle+\frac{s^{2}}{2 !}\left\langle x^{2}\right\rangle+\frac{s^{3}}{3 !}\left\langle x^{3}\right\rangle+\ldots \\[4pt] &=\sum_{n=0}^{\infty} \frac{s^{n}\left\langle x^{n}\right\rangle}{n !} \end{aligned} \nonumber \]

Para ordenar\(1 / N_{0}\), tenemos

\[\Psi(s)=\sum_{n=0}^{\infty} \frac{s^{n}\left\langle x_{0}^{n}\right\rangle}{n !}+\frac{1}{N_{0}} \sum_{n=0}^{\infty} \frac{s^{n}\left\langle x_{1}^{n}\right\rangle}{n !}+\mathrm{O}\left(1 / N_{0}^{2}\right) \nonumber \]

Ahora, usando (3.2.24),

\[\begin{aligned} \sum_{n=0}^{\infty} \frac{s^{n}\left\langle x_{0}^{n}\right\rangle}{n !} &=\sum_{n=0}^{\infty} \frac{\left(s e^{b t}\right)^{n}}{n !} \\[4pt] &=e^{s e^{b t}} \end{aligned} \nonumber \]

y usando\((3.2.27)\),

\[\begin{aligned} \sum_{n=0}^{\infty} \frac{s^{n}\left\langle x_{1}^{n}\right\rangle}{n !} &=\frac{1}{2}\left(1-e^{-b t}\right) \sum_{n=0}^{\infty} \frac{1}{n !} n(n-1)\left(s e^{b t}\right)^{n} \\[4pt] &=\frac{1}{2}\left(1-e^{-b t}\right) s^{2} e^{2 b t} \sum_{n=0}^{\infty} \frac{1}{n !} n(n-1)\left(s e^{b t}\right)^{n-2} \\[4pt] &=\frac{1}{2}\left(1-e^{-b t}\right) s^{2} \sum_{n=0}^{\infty} \frac{1}{n !} \frac{\partial^{2}}{\partial s^{2}}\left(s e^{b t}\right)^{n} \\[4pt] &=\frac{1}{2}\left(1-e^{-b t}\right) s^{2} \frac{\partial^{2}}{\partial s^{2}} \sum_{n=0}^{\infty} \frac{\left(s e^{b t}\right)^{n}}{n !} \\[4pt] &=\frac{1}{2}\left(1-e^{-b t}\right) s^{2} \frac{\partial^{2}}{\partial s^{2}}\left(e^{s e^{b t}}\right) \\[4pt] &=\frac{1}{2}\left(1-e^{-b t}\right) s^{2} e^{2 b t} e^{s e^{b t}} \end{aligned} \nonumber \]

Por lo tanto,

\[\Psi(s)=e^{s e^{b t}}\left(1+\frac{1}{2 N_{0}}\left(1-e^{-b t}\right) s^{2} e^{2 b t}+\ldots\right) \nonumber \]

Podemos reconocer el término entre paréntesis de (3.2.29) como una expansión de la serie Taylor de una función exponencial truncada a primer orden, es decir,

\[\exp \left(\frac{1}{2 N_{0}}\left(1-e^{-b t}\right) s^{2} e^{2 b t}\right)=1+\frac{1}{2 N_{0}}\left(1-e^{-b t}\right) s^{2} e^{2 b t}+\mathrm{O}\left(1 / N_{0}^{2}\right) \nonumber \]

Por lo tanto, a primer orden en\(1 / N_{0}\), tenemos

\[\Psi(s)=\exp \left(s e^{b t}+\frac{1}{2 N_{0}}\left(1-e^{-b t}\right) s^{2} e^{2 b t}\right)+\mathrm{O}\left(1 / N_{0}^{2}\right) \nonumber \]

Los libros estándar sobre teoría de probabilidad (e.g., Un primer curso de probabilidad de Sheldon Ross, página 365) detallan la derivación de la función generadora de momento de una variable aleatoria normal:

\[\Psi(s)=\exp \left(\langle x\rangle s+\frac{1}{2} \sigma^{2} s^{2}\right), \quad \text { for a normal random variable; } \nonumber \]

y comparar (3.2.32) con (3.2.31) nos muestra que la distribución de probabilidad\(P(x, t)\) a primer orden in\(1 / N_{0}\) es normal con la media y varianza dadas por

\[\langle x\rangle=e^{b t}, \quad \sigma_{x}^{2}=\frac{1}{N_{0}} e^{2 b t}\left(1-e^{-b t}\right) \nonumber \]

La media y varianza de\(x=N / N_{0}\) es equivalente a las derivadas para\(N\) en (3.1.25), pero ahora aprendemos que\(N\) es aproximadamente normal para poblaciones grandes.

La aparición de una distribución de probabilidad normal (también llamada distribución de probabilidad gaussiana) en una expansión de primer orden es de hecho un caso particular del Teorema del Límite Central, uno de los teoremas más importantes y útiles en probabilidad y estadística. Aquí señalamos una versión simple de este teorema sin pruebas:

Teorema de Límite Central: Supongamos que\(X_{1}, X_{2}, \ldots, X_{n}\) son variables aleatorias independientes e idénticamente distribuidas (iid) con media\(\langle X\rangle\) y varianza\(\sigma_{X}^{2}\). Entonces, para suficientemente grande\(n\), la distribución de probabilidad del promedio de los\(X_{i}\)'s, denotada como la variable aleatoria\(Z=\frac{1}{n} \sum_{i=1}^{n} X_{i}\), es bien aproximada por un gaussiano con media\(\langle X\rangle\) y varianza\(\sigma^{2}=\) \(\sigma_{X}^{2} / n\)

El teorema del límite central se puede aplicar directamente a nuestro problema. Considera que nuestra población está formada por\(N_{0}\) fundadores. Si\(m_{i}(t)\) denota el número de individuos descendientes del fundador\(i\) en el momento\(t\) (incluido el fundador aún vivo), entonces el número total de individuos en el momento\(t\) es\(N(t)=\sum_{i=1}^{N_{0}} m_{i}(t)\); y el número promedio de descendientes de un solo fundador es\(x(t)=N(t) / N_{0}\). Si el número medio de descendientes de un solo fundador es\(\langle m\rangle\), con varianza\(\sigma_{m}^{2}\), entonces aplicando el teorema del límite central para grandes\(N_{0}\), la función de distribución de probabilidad de\(x\) es bien aproximada por un gaussiano con media\(\langle x\rangle=\langle m\rangle\) y varianza\(\sigma_{x}^{2}=\sigma_{m}^{2} / N_{0}\). Comparando con nuestros resultados\((3.2.33)\), encontramos\(\langle m\rangle=e^{b t}\) y\(\sigma_{m}^{2}=e^{2 b t}\left(1-e^{-b t}\right)\).