Много написано математических моделей для экосистем "низшие растения - животные", "высшие растения - животные", но нет ни одной адаптированной к нашим задачам, а именно к системам жизнеобеспечения рыбы (очистка воды от продуктов жизнедеятельности рыб, снабжение рыб кислородом, утилизация углекислого газа), например, культивирования спирулины и хлореллы в фотореакторе. Конечно можно строить фотореактор, тратить деньги, а дальше ставить опыты и думать что получилось и почему. Но мы делаем всегда не так - сначала все рассчитываем, пытаемся спрогнозировать, смоделировать, а потом строим и убеждаемся, что все правильно рассчитали (в пределах допустимых ошибок).

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

Производительность биологической системы обеспечения жизнедеятельности осетров связана со структурой эксплуатируемой популяции организмов. Для решения различного рода задач необходимы математические модели, учитывающие отмеченную зависимость. Целью настоящей работы и являлось построение такой модели для системы жизнеобеспечения рыб, основанной на культивировании одноклеточной водоросли хлореллы.

При решении поставленной задачи мы учитывали неоднородность популяции только с точки зрения различной степени зрелости образующих ее клеток хлореллы. Согласно общепринятым представлениям, выделим четыре группы клеток, находящихся в первой, второй, третьей и четвертой фазах развития. Обозначим численности этих групп через x₁, x₂, x₃, x₄. Предположим, что все фазы имеют одинаковую длительность, равную dt. При условии, что все клетки благополучно заканчивают свой цикл развития и каждая из них производит в среднем z новых клеток, численности выделенных групп клеток в моменты времени tk и tk+dt связаны следующими уравнениями:

(1)

Система уравнений (1) является математической моделью идеальной популяции клеток хлореллы. В действительности длительность жизненного цикла и отдельных его фаз - случайные величины, поскольку отдельные клетки развиваются в не совсем одинаковых условиях и потенциальные возможности самих клеток различны. Поэтому за время dt не все клетки выделенной группы перейдут в следующую фазу развития. Часть клеток останется в той же фазе, а, следовательно, и в той же группе. Другая часть обгонит основную массу клеток и за время dt пройдет две фазы развития. Так, некоторые клетки первой группы на следующем временном шаге окажутся в третьей группе, а клетки второй группы - в четвертой. Быстро развивающиеся клетки третьей группы за время dt успеют пройти четвертую фазу и разделиться, дав определенное количество клеток в первую группу. А потомки некоторых клеток четвертой группы успеют за указанное время достигнуть второй фазы развития и оказаться во второй группе. Кроме того, за время dt часть клеток каждой группы погибает. С учетом различного темпа развития клеток и частичной их гибели система уравнений динамики численности популяции принимает следующий вид:

(2)

Через b₁₁ , b₂₁ и b₃₁ обозначено относительное количество клеток первой группы с замедленным, нормальным и ускоренным темпом развития соответственно. Аналогичный смысл имеют коэффициенты b₂₂, b₃₂, b₄₂ для второй, коэффициенты b₃₃, b₄₃ и b₁₃ для третьей и коэффициенты b₄₄, b₁₄ и b₂₄ для четвертой группы клеток.

Коэффициенты

(3)

характеризуют выживаемость клеток выделенных групп. Их значения лежат в пределах 0-1. Первое крайнее значение соответствует полной гибели клетки группы, а второе - отсутствию гибели.

Через Z₁₃ в системе (2) обозначено среднее число автоспор в одной клетке третьей группы с ускоренным темпом развития. Аналогичный смысл имеют коэффициенты Z₁₄ и Z₂₄ для клеток четвертой группы с нормальным и ускоренным темпом развития.

Запишем систему уравнений (2) в векторно-матричных обозначениях

где х(к) и х(к+1) - векторы размерности (4 * 1), характеризующие численность выделенных групп клеток в к-й и (к+1)-й моменты времени;

(5)

- матрица, элементы которой выражаются через коэффициенты уравнений системы (2) следующим образом:

(6)

Далее рассмотрим процесс отбора продукции эксплуатируемой популяции клеток хлореллы в системе.

Предположим, что слив части культуры из реактора и добавка питательной среды с указанной целью осуществляется импульсами, длительностью которых мы пренебрегаем. Допустим также, что разбавление культуры осуществляется с периодом, равным временному шагу нашей модели. Это означает, что через интервалы времени dt происходит мгновенное уменьшение численности популяции. При равномерном распределении клеток по объему реактора справедливо уравнение

где х(к) и у(к) - векторы численности выделенных групп клеток в к-й момент времени до и после разбавления культуры, а m - степень разбавления (скалярная величина).

Необходимое условие получения непрерывного процесса культивирования водорослей - стабилизация среды в реакторе, в том числе освещенности клеток. Последний параметр при прочих равных условиях определяется численностью самой популяции и ее структурой, так как клетки хлореллы в равных фазах развития обладают различными оптическими свойствами и размерами. Практически стабилизируется не освещенность, а оптическая плотность суспензии в реакторе, величина которой Р связана со структурой популяции векторным уравнением

где B - вектор размерности (1 * 4) коэффициентов, характеризующих влияние клеток выделенных групп на оптическую плотность культуры; х - вектор численности групп клеток.

Процессы отбора продукции популяции и регулирования освещенности в реакторе в реальных системах технологически связаны, т.е. в результате каждого разбавления культуры оптическая плотность снижается до заданной величины Р₃.

Таким образом, с учетом принятых допущений, динамику численности популяции клеток хлореллы и ее фазовой структуры можно описать системой векторно-матричных уравнений

(9)

где х и у - векторы численности выделенных групп клеток в к-й и (k+1)-й моменты времени до и после разбавления культуры соответственно.

Первое уравнение описывает изменение структуры популяции в результате процессов рождения, развития и гибели клеток между двумя разбавлениями (за время dt). Второе уравнение учитывает мгновенное уменьшение численности популяции при разбавлении культуры. Третье уравнение определяет степень разбавления культуры на каждом шаге (в конце интервала dt).

Появление одной клетки хлореллы в i-й группе сопровождается выделением в окружающую среду и поглощением из нее за время dt определенного количества веществ. Введем вектор g размерности         (1 * 4), характеризующий в указанном смысле массообмен выделенных групп клеток по одному веществу. Тогда общее количество этого вещества, поглощенное (выделенное) популяцией, можно определить векторным уравнением

Уравнения такого типа характеризуют метаболизм популяции. В частности, по ним может быть определена производительность системы по кислороду и углекислому газу между k-м и (k+1)-м моментами времени (за время dt).

Производительность системы по биомассе на этом же отрезке времени определяется векторным уравнением

где C - вектор размерности (1 * 4), характеризующий массу клеток выделенных групп, и m(k+1) - степень разбавления культуры в (k+1)-й момент времени.

Итак, мы получили математическую модель динамики производительности системы культивирования одноклеточной водоросли хлореллы. В стационарных режимах работы системы элементы матрицы А и компоненты векторов g, B и C - константы, зависящие от условий культивирования хлореллы. В нестационарных условиях указанные параметры системы являются также функциями времени.

В заключение необходимо затронуть вопрос практического определения численных значений параметров предложенной математической модели системы. Прежде всего, для этих целей могут быть использованы литературные данные по онтогенезу клетки хлореллы, содержащие сведения о длительности жизненного цикла клеток различных штаммов, количестве автоспор в одной материнской клетке, закономерностях накопления биомассы клетки и ее фотосинтеза по времени при разных условиях. Эти данные позволяют оценить численные значения временного шага модели dt и коэффициентов матрицы А, векторов g и С.

На основании такого рода данных нами была получена математическая модель производительности системы культивирования водоросли при выполнении ею функции регенерации кислорода в двух режимах работы - стационарного и режима, связанного с временным выключением источника света.

Были подобраны коэффициенты для матрицы А.

Она соответствует случаю, когда клетки не гибнут, каждая клетка делится в среднем на четыре автоспоры, а клетки с замедленным, нормальным и ускоренным темпом развития соотносятся, как 0,25 : 0,5 : 0,25.

Вектор g, в данном случае характеризует газообмен выделенных групп клеток в относительных единицах

Для процессов, происходящих в популяции при изменении светового режима, приняты следующие допущения: при переходе в темноту развитие клеток первой и второй групп приостанавливается, клетки третьей и четвертой групп продолжают развиваться, но делятся они на меньшее число автоспор, через время 2 dt развитие полностью прекращается, на каждом временном шаге (за время dt ) гибнет 10% клеток, при включении источника света клетки сразу же начинают развиваться, но на двух первых шагах делятся на меньшее число автоспор, гибель клеток прекращается.

Поэтому при выключении света происходит уменьшение скорости роста, остановка и гибель части водоросли. Тоже самое происходит при включении, но с другой динамикой. Поэтому были получены матрицы для уравнения (4) для темного периода и для светлого. Эти матрицы имеют зависимость во времени. Для того, чтобы не перегружать наших читателей мы опустим их описание. В приведенной ниже компьютерной программе Краснобородько В.В. используются эти матрицы для темного и светлого периода.

Вектор, характеризующий газообмен популяции в темноте g_Т , получен при условии, что интенсивность дыхания клеток составляет 5% интенсивности фотосинтеза, т.е. g_t= -0,05g.

Описанную модель использовали для теоретического анализа аварийных ситуаций в системе жизнеобеспечения, связанных с выключением на равное время источника света. На диаграммах показаны переходные процессы в производительности системы и численности клеток разных групп.

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

Существуют более сложные математические модели, их объяснение займет много времени и потребует специальных знаний высшей математики. Поэтому мы ограничимся описанием этой простой модели.

Написал программу и подготовил визуализацию Василий Краснобородько.

ОПРЕДЕЛЕНИЕ ПАРАМЕТРОВ МАТЕМАТИЧЕСКОЙ МОДЕЛИ ПОПУЛЯЦИИ ОДНОКЛЕТОЧНЫХ ВОДОРОСЛЕЙ ПРИ СТАЦИОНАРНЫХ УСЛОВИЯХ

На предыдущей странице мы рассмотрели математическую модель. Теперь экспериментальным путем найдем значения констант в этой модели, а также рассмотрим вопросы экспериментальной проверки эффективности методики идентификации параметров математической модели динамики размерной структуры популяции одноклеточных водорослей, совершенствования этой методики и оценки адекватности самой модели.

Математическую модель динамики размерной структуры популяции одноклеточных водорослей при стационарных условиях представим в следующем виде:

где x(k) и x(k+1) - векторы (n *1)- размерной структуры популяции в k-й и (k +1)-й момент времени; F - матрица (n* n)-популяции.

Идентификации (определению) подлежат элементы матрицы популяции F. Исходной информацией является регистрируемый процесс динамики структуры популяции одноклеточных водорослей при стационарных условиях. Для экспериментальных исследований была взята водоросль Chlorella vulgaris. Водоросль выращивали на питательной среде Тамийя. Опыты вели в экстенсивных условиях с алгологически чистой культурой в самодельном плоском культиваторе, описанном на этом сайте, объемом 40 литров. Исследования проводили в люминостате, где размещались колбы с суспензией, непрерывно барботируемой воздухом с содержанием С0₂ 5% . Температуру воздуха в люминостате поддерживали в пределах 37 ± 1°С. Освещенность равнялась 5000 лк. Опыт продолжался 60 ч, при этом начальная плотность популяции составила 4*I0³ а конечная 16*10⁶ кл/л.

Для измерения структуры водоросли исследуемую суспензию тщательно перемешивали и из разных точек объема колбы брали пробы, которые образовывали контролируемую суспензию объемом 10 мл. Эту контролируемую суспензию в свою очередь перемешивали и часть ее помещали в камеру Горяева.

Структуру водоросли измеряли путем фотографирования под микроскопом суспензии, находящейся в камере Горяева. При 40-кратном увеличении объектива микроскопа и 7-кратном объектива фотоаппарата, т.е. при общем 280-кратном увеличении при фотографировании, получена одинаковая и достаточно высокая резкость клеток по всей площади кадра. Фотографирование для повышения резкости граней клеток вели в поляризованных лучах фазово-контрастным методом.

Для определения структуры исследуемой водоросли клетки разбивали на пять размерных групп и подсчитывали число клеток в каждой группе. К первой группе относили клетки размером до 1,75 мкм, ко второй - от 1,75 до 2,25, к третьей - от 2,25 до 2,75, к четвертой - от 2,75 до 3,25 и к пятой группе - свыше 3,25 мкм. Распределение клеток по группам характеризует размерную структуру водоросли.

При такой методике (микрофотографирование) относительная погрешность измерения численности каждой группы структуры водоросли Chlorella vulgaris по результатам обработки одного фотокадра не превышает 90%. Увеличивая количество обрабатываемых фотокадров N, можно существенно уменьшить эту погрешность. Так, при N = 20 относительная погрешность измерения численности каждой группы структуры исследуемой водоросли не превышает  20%.

Эффективность методики идентификации зависит от начального состояния структуры исследуемой водоросли (равновесная или неравновесная структура) и шага измерения векторов структуры популяции. Если начальная структура популяции равновесная, то векторы структуры, измеренные в последовательные моменты времени с любым шагом измерений, близки к линейной зависимости и эффективность методики идентификации низкая. Если же начальная структура неравновесная, то первые измеренные векторы структуры популяции, при соответствующем выборе шага измерений, линейно независимы и для этих измерений методика эффективна.

Обратите внимание на графики роста популяции хлореллы на предыдущей странице. При начальных условиях концентрация клеток разных групп разная. При дальнейшем выращивании мелких групп становиться всегда больше крупных.

Не равновестность начальной структуры (вектор x(0)) и шага измерений dt =t (k+1)- t(k), (k=0,1,...) устанавливается критерием эффективности идентификации. На основании этого критерия можно утверждать, что методика идентификации эффективна, если среди коэффициентов

где xi- (0) и xi (1) - число клеток в i -й группе структуры водоросли, в первые два момента измерений найдется хотя бы один отрицательный и один положительный, абсолютные значения которых превышают относительную погрешность измерений векторов структуры. При выполнении этого критерия начальная структура водоросли будет неравновесной. Шаг измерений dt должен подбираться экспериментально.

Для получения неравновесной структуры водоросль, предварительно подращиваемая в люминостате, была перенесена в темноту. Измерять структуру по описанной выше методике начали спустя 12 ч после помещения водоросли в люминостат из темноты для того, чтобы произошла адаптация и стабилизировались параметры развития водоросли. Суспензию фотографировали каждые 2 ч (по 20 кадров, которые использовали для статистической обработки). Расчеты показали, что в наших опытах критерий (2) выполняется при dt= 4 ч. Динамика структуры хлореллы в эксперименте показана на рисунке.

Параметры математической модели динамики структуры водоросли (элементов матрицы популяции F) определяли методом регуляризации Тихонова. Эффективность этого метода может быть повышена применением статистического моделирования. Для этого сформируем исходное уравнение идентификации

В этом уравнении А - матрица измерений

(4)

f_i - вектор искомых элементов i-й строки матрицы F, а U_i - вектор измерений

и

(5)

Затем элементы матрицы и векторы измерений случайным образом изменялись внутри своих доверительных интервалов и для каждой реализации {A_l,U_il}, (l= ), где М - число испытаний, находили регуляризованное решение уравнения (3). Искомое решение определяли по формуле

(6)

Найденная с помощью рассмотренного способа матрица популяции имеет вид

(7)

Моделирование на компе показало, что при погрешности измерений, имеющей место в данном эксперименте, ошибка определения элементов матрицы F не превышает 60%. Для повышения точности определения элементов матрицы F были применены методы теории распознавания образов. Представим уравнение (1) в следующем виде:

Из этого уравнения следует, что вектор структуры популяции в k-й момент времени x(k) содержит информацию о векторе начального состояния x(0) и элементах матрицы популяции   F. Можно утверждать, что при одном и том же векторе начального состояния различным матрицам популяции в k-й момент времени соответствуют различные векторы структуры популяции.

Точное значение вектора начального состояния не известно, что связано о наличием ошибок измерений векторов структуры популяции. По экспериментальной информации может быть восстановлена функция распределения плотности вероятностей f(x(0)) этого вектора, подчиняющаяся нормальному закону распределения.

Будем считать, что конечное множество матриц f={F₁,F₂,..,F_R} достаточно полно характеризует область возможных значений элементов искомой матрицы популяции F. Эта область ограничена доверительными интервалами элементов матрицы (7). Тогда под образом матрицы Fj принадлежит f(j = i, R )   в пространстве векторов структуры популяции будем понимать совокупность функций распределения плотности вероятностей

также подчиняющихся нормальному закону распределения.

Повысить компактность образов можно путем кодирования функции распределения плотности вероятностей f(x_j(k)(j=i,R; k = 1,2,...). Для этого закодируем векторы структуры x (к), (k = 1,2,...) разложением по системе наиболее приспособленных ортонормированных векторов. Представим m-ю случайную реализацию вектора структуры популяции x_j^m(к) в следующем виде:

где Ψ_i(к)- ортонормированные векторы; С_ij^m(к)- коэффициенты разложения; r - количество членов разложения (r≤n). Наиболее приспособленные ортонормированные векторы являются собственными векторами матрицы

где М - математическое ожидание.

В результате кодирования каждой случайной реализации вектора X_j^m(k) размерности (nх1), принадлежащей j-му образу, будет соответствовать случайный вектор:

(12)

элементами которого являются коэффициенты разложения. Введем в рассмотрение расширенный вектор, который назовем кодовым

(13)

элементами которого являются векторы (12) в к последовательных моментах времени.

Тогда под образом матрицы   F_j (j =1,R)  в пространстве кодовых признаков будем понимать функцию распределения плотности вероятностей f(C_j)(j = 1,R). Размерность многомерной функции распределения f(C_j) выбирается экспериментально из условия получения устойчивого алгоритма классификации.

Задача классификации, т.е. собственно задача распознавания, рассматривается как задача проверкиR  статистических гипотез. Байесово решающее правило для функции распределения f(C_j) , подчиняющейся нормальному закону распределения, имеет вид

(14)

для всех j ≠i .

С помощью рассмотренных методов элементы матрицы F в условиях проведенного эксперимента могут быть определены с ошибкой, не превышающей 20%. Найденная матрица популяции имеет вид

(15)

В соответствии с этой матрицей, за дискретный интервал времени dt  = 4 ч в группах водоросли 85% клеток останется в первой группе, а 15% перейдет во вторую; 60% клеток останется во второй группе, 30% перейдет в третью, а 10% в четвертую группу; 55% клеток останется в третьей группе, 25% перейдет в четвертую и 5% в пятую, а 15% клеток разделится на две автоспоры, которые попадут в первую группу; 20% клеток останется в четвертой группе; 10% перейдет в пятую, а 70% клеток разделится в среднем на 2,21 автоспоры, причем 77% из них попадет в первую группу, а 23% - во вторую, в пятой группе останется 10% клеток, а 90% разделится в среднем на 3,89 автоспоры, причем 86% из них попадет в первую группу, а 19% - во вторую.

Адекватность найденной модели была оценена путем сравнения экспериментального и теоретического процессов динамики структуры водоросли на всем интервале наблюдения. Относительная среднеквадратичная ошибка отклонения сравниваемых процессов для каждой группы структуры водоросли не превысила 5%. Этот факт подтверждает адекватность построенной модели и эффективность предложенной методики идентификации.

По мат. http://www.krasnoborodko.com/test/chlorella_2.htm