Дифференциальные модели в экономике, биологии и медицине

В этом параграфе мы разберем несколько классических моделей, предложенных за последние 200 лет в различных областях науки: экономике, биологии и медицине. Общим при построении этих моделей является использование дифференциальных уравнений с разделяющимися переменными, которые несложно решить, прочитав §59 данного справочника.

п.1. Экономика. Равновесная цена в модели Вальраса

Начальные сведения о модели рыночного равновесия, кривых спроса и предложения – см. §18 справочника для 9 класса.
Рассмотрим поведение рыночной цены при небольшом отклонении от точки равновесия по методу, предложенному Леоном Вальрасом (1874 г.)
Пусть p - цена товара, D(p) - спрос на него, S(p) - предложение.
Пусть спрос и предложение на рынке уравновешены, равновесная цена равна \(p_0\).
Если спрос начнет немного превышать предложение, то цена начнет расти: $$ \frac{dp}{dt}=D(p)-S(p) $$ В точке равновесия: $$ \frac{dp}{dt}|_{p=p_0}=0,\ \ D(p_0)-S(p_0)=0 $$ Разложим каждую из функций с помощью дифференциала (см. §52 данного справочника) с точностью до линейного множителя: \begin{gather*} D(p)\approx D(p_0)+D'(p_0)(p-p_0)\\ S(p)\approx S(p_0)+S'(p_0)(p-p_0) \end{gather*} Тогда разность спроса и предложения: \begin{gather*} D(p)-S(p)\approx D(p_0)+D'(p_0)(p-p_0)-S(p_0)-S'(p_0)(p-p_0)=\\ =\underbrace{D(p_0)-S(p_0)}_{=0}+\left(D'(p_0)-S'(p_0)\right)(p-p_0)=\left(D'(p_0)-S'(p_0)\right)(p-p_0) \end{gather*} Получаем уравнение с разделяющимися переменными: \begin{gather*} \frac{dp}{dt}=\left(D'(p_0)-S'(p_0)\right)(p-p_0)\\ \frac{dp}{p-p_0}=\left(D'(p_0)-S'(p_0)\right)dt \end{gather*} Интегрируем: $$ \int \frac{dp}{p-p_0}=\ln(p-p_0),\ \ \int\left(D'(p_0)-S'(p_0)\right)dt=\left(D'(p_0)-S'(p_0)\right)t $$ Общее решение: $$ \ln(p-p_0)=\left(D'(p_0)-S'(p_0)\right)t+C $$ Пусть при \(t=0\) наблюдается неравновесная цена \(p(0)\ne p_0\). Находим C: $$ \ln(p(0)-p_0)=\left(D'(p_0)-S'(p_0)\right)\cdot 0+C\Rightarrow C=\ln(p(0)-p_0) $$ Решение задачи Коши: \begin{gather*} \ln(p-p_0)=\left(D'(p_0)-S'(p_0)\right)t+\ln(p(0)-p_0)\\ e^{\ln(p-p_0)}=e^{\left(D'(p_0)-S'(p_0)\right)t+\ln(p(0)-p_0)}\\ p-p_0=(p(0)-p_0)e^{\left(D'(p_0)-S'(p_0)\right)t} \end{gather*}

Например:
Пусть \(D(p)=9-\frac{p^2}{6},\ S(p)=\frac{p^2}{12}\)
\(\sqrt{y+1}=y'x\) - ДУ первого порядка первой степени
Тогда равновесная цена \(9-\frac{p_0^2}{6}=\frac{p_0^2}{12}\Rightarrow \frac{p_0^2}{4}=9\Rightarrow p_0^2=36\Rightarrow p_0=6\)

Значения производных: \begin{gather*} D'(p)=0-\frac{2p}{6}=-\frac p3,\ \ D'(p_0)=-2\\ S'(p)=\frac{2p}{12},\ \ S'(p_0)=1 \end{gather*} Изменение цены со временем в этом случае: $$ p(t)=6+(p(0)-6)e^{(-2+1)t}=6+(p(0)-6)e^{-t} $$ Построим графики для трех различных цен в начальный момент времени: $$ p(0)=\left\{5;7;9\right\} $$ \(p(0)=5:\ p(t)=6-e^{-t}\)
\(p(0)=7:\ p(t)=6+e^{-t}\)
\(p(0)=9:\ p(t)=6+3e^{-t}\)

Все три кривые постепенно сходятся к равновесной цене \(p_0=6\).

Устойчивое схождение к \(p_0\) будет наблюдаться только при условии: $$ D'(p_0)-S'(p_0)\lt 0 $$ Т.е кривая спроса должна быть более крутой в своем спуске, чем кривая предложения на подъеме. Говорят, что эластичность спроса по цене в точке равновесия должна быть выше, чем эластичность предложения по цене.

Если степень при экспоненте будет положительной, \(D'(p_0)-S'(p_0)\gt 0\) решение уходит на бесконечность. Говорят, что такое решение неустойчиво.
Если степень при экспоненте будет равна нулю, \(D'(p_0)-S'(p_0)=0\), цена не будет меняться и останется неравновесной.

п.2. Биология. Логистическое уравнение Ферхюльста для роста популяции

Пусть \(P(t)\) – численность популяции. Построим модель её изменения со временем.
Логично предположить, что прирост потомства в популяции пропорционален количеству особей, из чего получаем:

Решением этого уравнения будет \(P(t)=P_0e^{rt}\) - уходящая в бесконечность экспонента. Что заставило Мальтуса заявить о грядущем перенаселении планеты и потребовать жестких мер по ограничению рождаемости.
Впрочем, в 1804 г. население Земли достигло первого миллиарда, а сегодня, несмотря на многие неприятности, на планете живет уже в 7 раз больше.
Введем в уравнение некий уменьшающий рост популяции фактор – естественный или искусственный – пропорциональный квадрату \(P\).
Ферхюльст (1838 г.) предложил такую форму записи: $$ \frac{dP}{dt}=rP\left(1-\frac PK\right) $$

Это уравнение с разделяющимися переменными. Решаем: $$ \frac{dP}{P\left(1-\frac PK\right)}=rdt $$ Разложим дробь \(\frac{1}{P\left(1-\frac PK\right)}\) на сумму простейших дробей (см. §7 справочника для 8 класса). $$ \frac{1}{P\left(1-\frac PK\right)}=\frac{K}{P(K-P)}=\frac1P+\frac{1}{K-P} $$ Получаем: $$ \left(\frac1P+\frac{1}{K-P}\right)dP=rdt $$ Интегрируем: \begin{gather*} \int\left(\frac1P+\frac{1}{K-P}\right)dP=\ln P-\ln(K-P)=\ln\frac{P}{K-P}\\ \int rdt=rt \end{gather*} Общее решение: $$ \ln\frac{P}{K-P}=rt+C $$ Пусть начальная популяция при \(t=0\) равна \(P_0\). Находим C: $$ \ln\frac{P_0}{K-P_0}=r\cdot 0+C\Rightarrow C=\ln\frac{P_0}{K-P_0} $$ Решение задачи Коши: \begin{gather*} \ln\frac{P}{K-P}=rt+\ln\frac{P_0}{K-P_0}\\ \ln\frac{P}{K-P}-\ln\frac{P_0}{K-P_0}=rt\\ \ln\left(\frac{P}{K-P}\cdot\frac{K-P_0}{P_0}\right)=rt\\ \frac{P}{K-P}\cdot\frac{K-P_0}{P_0}=e^{rt}\\ P=e^{rt}\cdot\frac{P_0}{K-P_0}(K-P)\\ P\left(1+e^{rt}\cdot\frac{P_0}{K-P_0}\right)=e^{rt}\cdot\frac{P_0}{K-P_0}K\\ P=\frac{e^{rt}\cdot\frac{P_0}{K-P_0}K}{1+e^{rt}\cdot\frac{P_0}{K-P_0}}=\frac{KP_0e^{rt}}{K-P_0+P_0e^{rt}}=\frac{KP_0e^{rt}}{K+P_0(e^{rt}-1)} \end{gather*}

Например:
Пусть популяция растет со скоростью \(r=0,1\) тыс/год
Начальное количество особей \(P_0=1\) тыс
Максимальное количество, которое способна прокормить данная территория, \(K=10\) тыс

Модель показывает, что через 70 лет популяция займет всю нишу, и её рост фактически прекратится.
На начальном этапе преобладает r-стратегия: бурное размножение и короткая продолжительность жизни.
Исчерпание ресурсов заставляет переходить на K-стратегию: низкий темп размножения и долгую жизнь.

Экспериментально рост популяции по кривой Ферхюльста был подтвержден в лабораторных условиях для мух-дрозофил. В естественных условиях для животных – и тем более, в рамках социума для людей – закономерность нарушается.

п.3. Медицина. Модель развития эпидемии SIR

Традиционной моделью, описывающей процесс развития эпидемии, является модель SIR (Susceptible/Infected/Recovered), предложенная У. Кермаком и А. Маккендриком в 1927 г.
Вся популяция в модели делится на три группы:

• \(S(t)\)— восприимчивые к инфекции, здоровые на момент времени \(t\);
• \(I(t)\)— уже инфицированные;
• \(R(t)\)— выздоровевшие, больше невосприимчивые к инфекции.

Популяция считается постоянной, т.е. \(N=S(t)+I(t)+R(t)=const\).
Рождаемость и смертность не учитывается.
Получаем следующую систему дифференциальных уравнений: $$ \begin{cases} \frac{dS(t)}{dt}=-\beta\frac{S(t)I(t)}{N}\\ \frac{dI(t)}{dt}=\beta\frac{S(t)I(t)}{N}=-\gamma I(t)\\ \frac{dR(t)}{dt}=\gamma I(t) \end{cases} $$ где \(\beta\) – скорость заражения, вероятность заболевания в случае контакта с инфицированным; \(\gamma\) - скорость выздоровления, \(\gamma=\frac1T;\ Τ\) – период болезни.
Начальные условия в момент времени \(t=0\): $$ S(0)=S_0\geq 0,\ \ I(0)=I_0\geq 0,\ \ R(0)=R_0\geq 0 $$ Переход из одной группы в другую можно изобразить линейной схемой:

Полученная система уравнений не является линейной и не имеет точного аналитического решения. Но её можно решить с использованием численных методов.
$$ \begin{cases} \frac{\triangle S(t)}{\triangle t}=-\beta\frac{S(t)I(t)}{N}\\ \frac{\triangle I(t)}{\triangle t}=\beta\frac{S(t)I(t)}{N}-\gamma I(t)\\ \frac{\triangle R(t)}{\triangle t}=\gamma I(t) \end{cases} $$ Считаем \(\triangle t=1\) – следующий шаг итерации. Тогда: $$ \begin{cases} S_{i+1}-S_i=-\beta\frac{S_iI_i}{N}\\ I_{i+1}-I_i=\beta\frac{S_{i+1}I_i}{N}-\gamma I_i\\ R_{i+1}-R_i=\gamma I_{i+1} \end{cases} $$ Получаем следующий итеративный процесс: $$ \begin{cases} S_{i+1}=\left(1-\beta\frac{I_i}{N}\right)S_i\\ I_{i+1}=\left(1+\beta\frac{S_{i+1}}{N}-\gamma \right)I_i\\ R_{i+1}=\gamma I_{i+1}+R_i \end{cases} $$ Знаний по информатике вам должно хватить, чтобы написать небольшой скрипт с циклом для этих уравнений и построить график.

Например:
Пусть общее количество населения N=10 тыс.чел.
В начальный момент инфицирован 1% населения: $$ S(0)=0,99N,\ \ I(0)=0,01N,\ \ R(0)=0 $$ Параметры: \(\beta=0,128;\ \gamma=0,096\) в расчете на день (эти параметры были рассчитаны по фактическим данным для лихорадки Эбола в Сьерра-Леоне).
Результат моделирования в MATLAB:
Красная кривая – это количество болеющих в данный момент. Как мы видим, к концу года она стремится к 0. Пик приходится на 70-80 дней с начала эпидемии и составляет 413 чел. или 4,13% населения.
Зеленая кривая – количество переболевших, к концу года выходит на асимптоту в 4700 чел. или 47,0% населения.
Синяя кривая – количество так и не заболевших, к концу года спускается на асимптоту в 5300 чел. или 53,0% населения.

Чем больше больных у вас будет в начале эпидемии и чем больше параметр \(\beta\), тем выше будет пик \(I_{max}\) для болеющих. Также, количество переболевших в конце эпидемии будет больше количества не заболевших.

Модель SIR – это начальный этап для исследований. На практике для моделирования эпидемий могут использоваться модели с десятками переходов и параметров, с постепенным усложнением по мере накопления данных.