Линеаризация информативных сигналов в микроаналитических приборах и методы их обработки

Размещено на http://www.allbest.ru/

Размещено на http://www.allbest.ru/

АВТОРЕФЕРАТ

диссертации на соискание ученой степени доктора физико-математических наук

ЛИНЕАРИЗАЦИЯ ИНФОРМАТИВНЫХ СИГНАЛОВ В МИКРОАНАЛИТИЧЕСКИХ ПРИБОРАХ И МЕТОДЫ ИХ ОБРАБОТКИ

Специальность 01.04.01 - «Приборы и методы экспериментальной физики»

Буляница Антон Леонидович

Санкт-Петербург - 2008

1. Актуальность темы

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

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

Реализуемое в последние годы направление миниатюризации приборостроения, в частности, сопровождается появлением исходных информативных сигналов новой структуры. Например, перевод ряда сепарационных методов химического и биологического анализа (электрофореза, хроматографии и т.п.) на микрочипы привел к изменению формы аналитических пиков от гауссовой к кусочно-линейной (трапециидальной или треугольной, т.е., совокупности ЛТ0 и ЛТ1).

Таким образом, основания для выбора в качестве унифицированного преобразованного сигнала совокупности ЛТ0 и ЛТ1: 1) это естественная форма исходного информативного сигнала ряда микроаналитических приборов; 2) малое число оцениваемых параметров: либо величина ЛТ0, либо, как правило, только параметр положения a ЛТ1 (b обычно связан с фоновыми измерениями); 3) простота обработки - интерполяция, дифференцирование и т.п.

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

Предпосылкой для унификации информативного сигнала к совокупности ЛТ1 может служить методика [1] идентификации 7-ми типов функциональных зависимостей. Этими зависимостями адекватно аппроксимируются информативные аналитические сигналы различных приборов. Тип зависимости идентифицируется на основе сравнения средних (арифметическое, геометрическое, гармоническое) входной и выходной величин, а зависимость приводится к линейному тренду с помощью преобразований, представленных в таблице.

Таблица 1 Функциональные зависимости и метод их преобразования к линейному тренду

Тип	Зависимость Y=F(X)	Метод перехода к виду Z=A+Bt	Алгоритм идентификации типа
1	Y=AX+B	Z=Y, t=X
2		Z=ln(Y), t=ln(X)
3		Z=ln(Y), t=X
4		Z=Y, t=1/X
5		Z=1/Y, t=X
6		Z=1/Y, t=1/X
7		Z=Y, t=ln(X)

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

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

Цель работы:

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

Достижение указанной цели обеспечивается решением следующих задач:

1. Исследовать структуру исходных (предварительных) информативных сигналов микроаналитических приборов, в т.ч. химического, биологического, иммунного анализа.

2. Обосновать форму унифицированного сигнала, как совокупности ЛТ0 и ЛТ1, и рассмотреть методы унификации к указанному виду с помощью простых нелинейных преобразований временной оси и предварительного информативного сигнала. Как следствие, определить класс микроаналитических приборов, допускающих подобную процедуру унификации.

3. Разработать метод оценивания параметров информативных сигналов типа ЛТ0 и ЛТ1, обладающий гарантированной эффективностью, робастностью, несмещенностью и состоятельностью оценок при экономичности и простоте реализации.

4. Исследовать применимость алгоритма стохастической аппроксимации Роббинса-Монро [2] и его модификаций [3-6] в качестве основы метода оценивания параметров унифицированного информативного сигнала.

5. Реализовать в форме программного продукта алгоритм оценивания параметров унифицированных информативных сигналов (параметра положении ЛТ1 или/и величины ЛТ0), также включающий выбор начальных установок (величины зоны нечувствительности, масштабного поправочного множителя, начального приближения оценки), и предусматривающий остановку оценивания при наличии разладки в последовательности измерений.

6. Решить ряд задач по предварительной обработке информативных сигналов, включающий: а) правило Новицкого-Зограф [7] для разработки унифицированного критерия отбраковки выбросов, б) оптимизация ширины медианного окна для повышения эффективности получения робастных оценок математического ожидания, в) анализ эффективности применения смещенных экстремальных порядковых статистик для оценивания математического ожидания в условиях аддитивной ограниченной помехи, г) цифровая фильтрация ЛТ1, для повышения эффективности оценивания параметров аналитических пиков (временное положение, амплитуда и площадь пика) в условиях малости отношения сигнал/шум.

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

Научная новизна

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

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

3. Подтверждены основополагающие свойства оценки предложенного алгоритма - несмещенность, состоятельность и квази-эффективность, понимаемая как решение минимаксной задачи (наименьшая дисперсия при наихудшем законе распределения помехи). Предложены новые методы доказательства несмещенности оценки: на основе модифицированного апостериорного критерия Аоки и с помощью интерпретации алгоритма нелинейной дискретной САУ, и анализа ее устойчивости.

4. Впервые выявлен единственный класс аддитивных помех (Симпсоновская помеха), для которого дисперсия оценки величины ЛТ0 при применении алгоритма Фабиана-Цыпкина не зависит от величины выбранной зоны нечувствительности.

5. Сформулированы требования к построению линейных калибровочных функций на основе концепции Хьюбера об исключении точек риска.

6. Проанализирован информативный сигнал прибора серии АНК, реализующего ПЦР в реальном масштабе времени в форме кривой логистического роста, доказана адекватность подобной формы кинетической кривой и показан метод сведения информативного сигнала к форме ЛТ1 с последующим оцениванием коэффициента увеличения биомассы пробы за 1 цикл.

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

8. Исследована эффективность алгоритмов первичной обработки информативных сигналов, использующих медианные порядковые статистики. Впервые решена задача оптимизации ширины медианного окна на основе минимизации показателей качества - финальной ошибки прогнозирования (ФОП) или/и информационного критерия Акаике (ИКА).

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

Практическая ценность работы.

1. Исследованы области применимости универсальных алгоритмов первичной обработки измерений: а) правило Новицкого-Зограф отбраковки выбросов, б) применение L-оценок на основе медианных порядковых статистик (выбор оптимального медианного окна); в) предложены новые алгоритмы оценивания площади электрофоретического пика при малости отношения сигнал/шум.

2. Модификация информативного сигнала к совокупности ЛТ1 позволяет унифицировать его обработку и базовое ПМО для широкого класса микроаналитических приборов

3. В качестве универсального алгоритма обработки сигналов типа ЛТ1 с оцениванием его параметра положения при аддитивной случайной помехе с априорно неизвестным законом распределения предложена комбинация перехода к первой разности с симметризацией помехи, и применение рекурсивного алгоритма стохастической аппроксимации. Алгоритм реализован в виде программного продукта.

4. Указанный алгоритм сигналов использован в базовом ПМО разработанных в Институте аналитического приборостроения РАН приборах: а) хемосенсорных анализаторах рН, ионов тяжелых металлов и редкоземельных элементов в водных средах (фотометры серий SEN и mSEN), б) прибора для фотоплетизмографического определения степени кислородного насыщения артериальной крови (CADIX OXI), в) приборах для определения числа и последовательности нуклеиновых кислот серий АНК 16 и АНК 32, г) микрофлюидных электрофоретических устройствах для анализа биологических проб.

5. Сформулированы требования построения линейных калибровочных функций, базирующиеся на необходимости исключения точек риска по Хьюберу.

Положения, выносимые на защиту.

1. Применимость кусочно-линейного сигнала с аддитивной случайной помехой в качестве унифицированного информативного сигнала прибора (принцип объединения приборов по принадлежности исходного сигнала к 7-ми типам и методы их трансформации к унифицированному виду).

2. Метод оценивания параметров кусочно-линейных сигналов (совокупности ЛТ1) при аддитивной помехе с априорно неизвестным законом распределения, и реализация соответствующего метода в форме программного продукта.

3. Доказательства несмещенности оценок методами: а) анализ устойчивости соответствующей нелинейной САУ (по критерию типа Попова), б) модифицированный статистический апостериорный анализ Аоки.

4. Доказательство существования единственного класса невырожденных аддитивных помех (треугольная Симпсоновская помеха), при наличии которых дисперсия ошибки оценивания ЛТ0 при применении алгоритма Цыпкина не зависит от величины зоны нечувствительности.

5. Критерий оптимизации выбора ширины медианного окна, что позволяет повысить эффективность первичной обработки информативных сигналов при сохранении робастности.

6. Модель отклика принципиально нового чувствительного элемента биосенсора - самоорганизующейся колонии несовершенных грибов и ее свойства: а) механизм самоорганизации и его управляющие параметры, б) стратегии развития колонии, в) способы формализации особенностей развития, г) интерпретация смен стратегий развития с позиций фазовых переходов 1-го и 2-го рода, д) представление информативных сигналов в форме совокупности ЛТ0 и ЛТ1.

7. Принцип построения линейных калибровочных функций, базирующийся на концепции Хьюбера об исключении точек риска.

Апробация результатов работы

Основные результаты работы докладывались и обсуждались на конференциях и научных школах: VI Межд. Конференция и Выставка по Инструментальному анализу (1995, Пекин, Китай), II и III Съезды биофизиков России (1999, Москва, 2004, Воронеж), 4-я Санкт-Петербургская Ассамблея молодых ученых и специалистов (1999), 2-й Всеросс. Симпозиум по Проточному Химическому анализу (1999, Москва), 2-я, 3-я и 4-я Межд. Конференции «Цифровая обработка сигналов и ее применения» (1999, 2000, 2002, Москва), Межд. конференция по комплексным системам (New England Complex Systems Institute - NECSI) (2000, Нашуа, США), Межд. семинар «Нелинейная динамика в биологии» (2000, Копенгаген, Дания), Росс. научно-практическая конференция «Оптика и научное приборостроение - 2000» ФЦП «Интеграция» (2000, Санкт-Петербург), 8-я и 9-я Межд. Конференции «Математика. Компьютер. Образование» (2001, Пущино, 2002, Дубна), LVI Научная сессия, посвященная Дню радио (2001, Москва), Int. Conference on Biological Physics (2001, Kioto, Japan), 1-й Съезд микологов России «Современная микология в России» (2002, Москва), Межд. конференция «Instrumentation in ecology and human safety» (2002, Санкт-Петербург), 6-я и 7-я Всеросс. школы-конференции «Биология - Наука 21 века» (2002, 2003, Пущино), 1-й и 2-й Межд. Симпозиумы «Биокосные взаимодействия: жизнь и камень» (2002, 2004, Санкт-Петербург), XXII Dynamics Days Europe (2002, Гейдельберг, Германия), 1-я Всеросс. конференция «Аналитические приборы» (2002, Санкт-Петербург), 2^nd Black Sea Conference on Analytical Chemistry, Workshop on 1-st Marmara Analytical Chemistry (2003, Стамбул, Турция), 3-й Межд. Симпозиум по сепарации в бионауках (2003, Москва), Межд. Форум «Аналитика и аналитики» (2003, Воронеж), 2-й и 3-й Всеросс. Конгрессы по медицинской микологии. (2004, 2005, Москва), Межд. Специал. Симпозиум ISSY25 (Systems biology of yeasts - from models to applications) (2006, Ханасаари, Финляндия), 8-й Всеросс. Симпозиум по прикладной и промышленной математике (2007, Адлер), на совместных Семинарах Института аналитического приборостроения и Химического факультета Московского государственного университета «Микрофлюидные технологии и их применения» (2005, 2007, Санкт-Петербург), на Семинарах кафедры Биофизики СПбГУ и Института аналитического приборостроения РАН.

Поддержка работы Программами и Грантами

Работы поддерживались Российскими и Международными Программами и Грантами, среди которых ФЦП «Интеграция» с участием СПб государственного университета (кафедра биофизики), Российско-Греческий проект «Развитие методов для документации и анализа данных по факторам биоразрушения музейных объектов», INTAS (1993 - 1996 гг.), Гранты поддержки молодых ученых Санкт-Петербурга (1998-1999 гг.), Индивидуальные Трэвел-гранты (2000, 2001 гг.), Программы СПб Научного центра РАН (с 2000 г.), МНТП «Вакцины нового поколения и диагностические системы будущего», Грант РФФИ № 03-01-39003ГФЕН_а «Теоретические и экспериментальные исследования явлений переноса и взаимодействия биологических объектов в микрофлюидных устройствах» и др.

Направления исследований тесно связаны с приоритетными направлениями развития науки, технологий и техники Российской Федерации (Информационно-телекоммуникационные системы; Живые системы) и критическими технологиями РФ (Биоинформационные технологии, технологии обработки, хранения, передачи и защиты информации, Технологии биоинженерии), утвержденными Президентом РФ 21.05.2006.

Личное участие автора

Вклад автора в полученные результаты диссертационной работы в областях теории обработки информации и методов построения калибровочных функций был определяющим. При разработке различных приборов (хемосенсорные анализаторы серии SEN, прибор определения кислородного насыщения крови CADIX, прибор для осуществления полимеразной цепной реакции (ПЦР) в реальном времени, микрофлюидная аналитическая система и др.) теоретические разработки с участием автора в области математического (имитационного) моделирования определили пути технического решения прибора, включая оптимизацию конструкции и режима измерения. На стадии практической реализации прибора вклад автора, главным образом, состоял в разработке теоретических основ ПМО, позволяющего производить эффективную обработку и оценивание параметров информативных сигналов. Автор возглавлял и реализовал направление теоретических исследований, связанное с разработкой концепции нового типа сенсорной системы (на основе самоорганизующейся колонии несовершенных грибов), включая создание математической модели базовых процессов, ее модификации, программной реализацией математического моделирования и разработкой алгоритмов оценивания. Лабораторные эксперименты проводились совместно с сотрудниками кафедры Биофизики СПбГУ (д.б.н. Л.К. Панина, к.б.н. Е.В. Богомолова, аспиранты Е.Ю. Быстрова, Е.О. Цветкова и И.А. Зароченцева).

Общее число публикаций - 66, включая 28 статей (из них, в журналах, внесенных в перечень ВАК РФ - 18), материалы более 30 конференций, 4 научных отчета.

2. Содержание работы

Во Введении обоснована актуальность темы, сформулированы цель и задачи, научная новизна, практическая ценность. Даны сведения об апробации работы, защищаемых положениях, основных публикациях и личном вкладе автора, описаны структура и объем диссертации.

Глава 1 посвящена анализу информативных сигналов различных аналитических приборов.

В первой части Главы исследуются экстракционные методы и приборы. Теоретические основы были заложены чл.-корр. Г.А. Ягодиным и его коллегами [8]. В диссертации для более детального обсуждения были выбраны 12 схем, для которых проведено обобщение результатов при произвольном значении параметров (константы скоростей реакций, коэффициенты диффузии и т.п.). Для каждой схемы была предложена математическая модель в форме системы дифференциальных уравнений реакция-диффузия. Эти уравнения решались методом Лапласа, а переход от изображений количества экстрагированного продукта Q(p) к оригиналам выполнен с помощью численного метода Папулиса, что позволило получить требуемые оценки при малых временах. Подтвердилось, что зависимость количества экстрагированного продукта q от времени t, начиная с некоторого момента, определяемого характеристиками реакций и схемой экстракции, адекватно описывается ЛТ1 в координатах , где a=1/2, реже 1. При этом, параметр положения ЛТ1 прямо пропорционален концентрации анализируемого вещества. Также рассматривалось формирование информативного сигнала в приборах химического экспресс-анализа с фотометрическим детектированием.

Во второй части описана процедура определения степени кислородного насыщения артериальной крови и частоты пульса на основе фотоплетизмографических измерений с помощью разработанного Институтом аналитического приборостроения РАН измерителя CADIX OXI. Отмечено, что важнейшей составляющей сигнала (постоянный оптический сигнал на 2-х длинах волн), является ЛТ0. Отслеживание этой составляющей (и ее компенсация) позволили повысить точность оценки переменной составляющей сигнала. Степень кислородного насыщение крови рассчитана как , где I_v1, I_v2 и I_c1, I_c2 - переменные (v) и постоянные (c) составляющие интенсивности отраженного от биоткани излучения; ₁ и ₂ - длины волн ИК и красного спектральных диапазонов (обычно, изобестические точки 650 и 805 нм); A и В - эмпирические коэффициенты. Поскольку амплитуда переменных составляющих сигнала многократно меньше постоянных, оценивание и компенсация последних позволяет измерить переменные составляющие с большой чувствительностью и точностью.

Еще одним подклассом аналитических приборов являются иммунные анализаторы, применяющие метод ПЦР в реальном масштабе времени для анализа количества фрагментов ДНК. Показано, что модель логистического роста Ферхюльста-Пирла адекватно описывает динамику информативного сигнала. На участке (t0) была проведена компенсация фонового сигнала, близкого по форме к ЛТ1. На стадии экспоненциального роста, удаленной от области насыщения, в полулогарифмических координатах сигнал линеаризован, и его параметр положения показывает коэффициент роста количества ДНК за 1 цикл ПЦР (при теоретическом пределе 2, он составил 1.78-1.83). С другой стороны, линеаризация информативного сигнала q(t) для всего временного диапазона выполняется преобразованием , где q_max - величина сигнала насыщения, априорно неизвестного.

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

Математическая модель конвективно-диффузионного массопереноса вещества на базе уравнения Навье-Стокса в микроканале чипа (совокупности плоских щелей), при маршевой координате x, координате сечения y и времени t, примет вид (1):

(1)

Здесь C=C(x,y,t) распределение концентрации, U(y) - профиль конвективной скорости. Граничные условия: симметрия на оси канала, непроницаемость стенки, недоступность вещества на бесконечности. Введение моментов порядка n: , где m₀ - количество вещества, h=m₁/m₀ - центр тяжести пика, дисперсия пика , позволило (1) свести к дифференциальному уравнению 2-го порядка для моментов вида (2):

(2).

Здесь z, u* и D* - безразмерные координата сечения, средняя конвективная скорость и коэффициент диффузии. Начальное условие: равномерно заполненная пробка вещества длины 2D, Pe=u*/D* -- число Пекле. Средние центр тяжести и дисперсия пика и , где B^j - коэффициенты разложения скоростного профиля () в ряд Фурье по косинусам. Третий и четвертый моменты найдены численно, что позволило оценить коэффициенты асимметрии и эксцесса пика, и подтвердить гипотезу о его кусочно-линейной форме. Главные тенденции изменения характеристик пика: а) стремление коэффициента асимметрии к нулю с увеличением времени анализа, б) коэффициент эксцесса стремится к 3, но при малых t он равен 2.3-2.6. Последнее более соответствует гипотезе треугольной или трапециидальной формам пика, а не гауссовой.

В целом, в Главе 1 рассмотрены некоторые типы информативных сигналов в аналитических приборах. Несмотря на различные принципы действия, приборы можно объединить в один класс, т.к. 1) их исходный информативный сигнал описывается одним из рассмотренных семи типов функциональных зависимостей, 2) все эти сигналы могут быть преобразованы к ЛТ1 или их совокупности с помощью трансформаций измеряемой величины или/и временной оси. Следовательно, информативный сигнал широкого класса микроаналитических приборов может быть преобразован к следующему виду: ЛТ1 (в специальной системе координат, определяемой типом прибора) с аддитивной случайной помехой с априорно неизвестным законом распределения. В этом случае информативным параметром является величина параметра положения.

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

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

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

Математическая модель самоорганизации в колониях мицелиальных грибов основывается на следующих положениях: 1. Колония является единым образованием, развивающимся радиально путем приращения и ветвления нитевидного мицелия по поверхности гелеподобной cреды, за счет потребления субстрата (активатора) и выделяющим продукты метаболизма (ингибитора роста). 2. Сценарии развития колоний: равномерный газон, кольцевые структуры (КС) с областями разреженности и сгущения, «умирание» колонии (локализация вблизи зоны начального посева) 3. Распределение субстрата и метаболитов происходит по механизму типа реакция-диффузия. 4. Коэффициент диффузии метаболитов значительно превышает коэффициент диффузии субстрата, т.е. предполагаем диффузию жидких метаболитов по условно твердой гелевой поверхности, и диффузию субстрата внутри пористой (условно твердой) гелевой среды.

Математическая модель пространственного роста колонии мицелиальных грибов имеет вид системы (3а)-(3в), описывающей распределение плотностей мицелия (), концентраций субстрата (s) и продуктов метаболизма (m):

. (3а)

(3б)

(3в)

Входящие в систему (3а)(3в) переменные и функции описаны в диссертации. Символ * в уравнении (3в) подразумевает «полную частную» производную , что позволило учесть локальный прирост плотности мицелия в заданной точке с долей и радиальное распространение нити с долей (1-). Учтен принцип инерционности: по достижении критической концентрации ингибиторов (m>⁰) с задержкой ⁰ происходит прекращение роста мицелиальной клетки. Другие особенности модели: 1. Граничные условия 2-го рода (непроницаемость и радиальная симметрия). 2. Скорость выработки продуктов метаболизма определяется функцией U₁. Рассмотрено 3 гипотезы: а) постоянство , б) старение , в) «запрет по насыщению с уровнем » , где - аналог предельного репродуктивного возраста. Сравнение гипотез (б и в) показало схожие радиальные распределения плотности мицелия (равное число максимумов при их пространственном смещении на 0,08-0,10 радиуса). Последнее объясняется тем, что в цикле роста мицелия максимумы выработки метаболитов для разных моделей смещены друг относительно друга.

Способность мицелиальных грибов реагировать на изменения в окружающей среде изменением формы колонии, в частности, стремление части мицелия «уйти» в области с более благоприятными условиями, учтена в модели по схеме: , где m - концентрация продуктов метаболизма, m⁰ - концентрация, вызывающая отравление. Т.о. при росте m доля мицелия (1-), покидающего зону отравления, будет увеличиваться. Учет этого эффекта адаптации существенен при переходном состоянии колонии между умиранием и КС, и он позволяет значимо сместить эту границу выживания в условия, когда отравление среды происходит быстрее.

Мерой упорядоченности колоний несовершенных грибов является энтропия, и измерение колонии на информационной основе связано с общими понятиями хаоса и порядка. Управляющий (по Климонтовичу) [9] параметр а влияет на структуру колонии следующим образом: его монотонное изменение приводит к монотонному изменению структуры от абсолютного хаоса при а=0 к идеальному порядку а=а_мах. Радиальное распределение плотности мицелия _i=(_i) в момент времени =_j после нормализации можно считать аналогом дискретной плотности распределения вероятностей и использовать для расчета энтропии. Аналогией понятия хаос является структура сплошной газон, наиболее упорядочена умирающая колония.

Сходство между экспериментально полученной структурой и модельным распределением плотности мицелия иллюстрируется рис. 1 (А и Б), что подтверждает адекватность модели (3а) - (3в).

А Б

Рис. 1. Фотография фрагмента колонии A.alternata с образованием КС (А) и соответствующее модельное распределение плотности мицелия (Б) при =0.90, v=1.2, =5, =0.4, s⁰=20, =0.2, D_m=1.10^{-3, 0=}0.1, =0.4, ⁰⁼0.50, ₀=0,03

Оценивание коэффициента диффузии метаболитов D_m позволило преобразовать модель (3а)-(3в) к размерным параметрам. За основу принято уравнение распределения метаболитов (3а) («моделью точечного источника»), т.к. основная часть метаболитов вырабатывается в области, где плотность мицелия наибольшая (изначально в точечной зоне посева в центре колонии).

Следствием радиальной диффузии будет приближенное представление концентрации продуктов метаболизма как

. (4)

При граничных условиях 2-го рода и начальном условии равномерного посева: и . Здесь - цилиндрические функции Бесселя нулевого и первого порядков, r - радиальная координата, t - время, соответствует стенке чашки Петри. Выбрав три первых слагаемых из (4), можно из отношения оценить и величину D_m. Получены оценки 0,200,55^.10^-6 см²/сек, что соответствует случаю диффузии низкомолекулярных соединений в жидкой среде.

Математическая модель кооперативного развития двух диморфных форм - необходимый этап, т.к. мицелий - не единственная клеточная форма. Модель диморфного (мицелий-дрожжи) перехода базировалась на модели Лотки-Вольтерра (аналог модели хищник-жертва). После нормировок концентраций клеточных форм (u_1,2) и времени () модель приведена к виду:

.(5)

Получены адекватные стационарные решения (5) - устойчивый узел и седловые точки. При различных соотношениях устойчивому решению соответствуют или полное доминирование одной из форм, или присутствие обеих клеточных форм. При (5) имеет два устойчивых состояния: их выбор определяется начальными концентрациями. Выявлена связь между параметрами модели (5) k_2,3 и внешними условиями (концентрациями марганца С_1, мг/л и глюкозы С_2, г/л) в форме неканонических эллиптических кривых:

При этом, модель (5) неадекватно описывает динамику диморфных переходов, допуская области «отмирания» . Была проведена ее модификация схожего вида, свободная от указанного недостатка:

.

В следующем разделе с позиций термодинамической (ТД) теории фазовых переходов Л.Д. Ландау рассмотрена самоорганизация мицелиальной формы колонии. Ранее [10] с этих позиций был исследован диморфный переход (переход от мицелиальной к дрожжевой форме роста).

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

переход стратегии от сплошного газона к КС может трактоваться с позиций фазового перехода 2-го рода. Роль ТД температуры играет - интенсивность роста массы мицелия. Ее критическое значение существенно меняется в зависимости D_m.

переход от «умирания» к КС не может быть описан в рамках двухфазной системы. Трактовка КС как метастабильного состояния (с различным расположением колец) позволяет описать последовательность переходов от одного состояния к другому с позиций фазовых переходов. Параметром порядка является энтропия, а переход от одного состояния к другому в рамках КС характеризуется ее скачками, что свойственно переходам 1-го рода.

В заключительном разделе Главы анализируется концепция применения сенсора на основе колонии несовершенных грибов для биотестирования окружающей среды. В частности, рассмотрена его связь с мультисенсорными системами «электронный нос» и «электронный язык». Колония несовершенных грибов интерпретирована как чувствительный элемент сенсора. Показано, что перспективно рассматривать колонию в целом как единый объект. Выявлено, что радиальные распределения плотностей мицелия при трех основных сценариях развития представляют собой зависимость, адекватно аппроксимируемую совокупностью ЛТ0 и ЛТ1.

Глава 3 посвящена анализу методов предварительной обработки информации. Проблема выявления грубых погрешностей рассматривается с позиций унифицикации алгоритма отбраковки выбросов. Достаточной степенью общности обладает утверждение Новицкого-Зограф - так называемое «замечательное свойство квантили 95%» - Z_0,95:

Z_0,95=m+k^. , (6),

где коэффициент k=1,650,05, m - математическое ожидание, ² - дисперсия случайной величины. Декларирована истинность (6) практически для всех наиболее распространенных законов распределения случайной величины. На плотность распределения вероятностей (ПРВ) накладываются ограничения: а) симметрия, б) высокий энтропийный коэффициент К_e>1,7, в) малая антимодальность (С_д1,5), где С_д=X_max/ при моде X_max. В этом случае правило отбраковки выбросов будет . Оно включает в отбрасываемые измерения не более 10% объема выборки (как полагаем, все измерения, соответствующие грубой погрешности). В диссертации проверена истинность (6) для наиболее распространенных одномодальных (около 20 случаев) и некоторых двумодальных распределений.

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

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

Альтернативный подход состоял в оценивании математического ожидания на основе смещенных экстремальных порядковых статистик, как это сделано в методе ПИО (простого интервального оценивания), получившего развитие в области анализа больших массивов физико-химических данных [11]. Одним из его этапов является интервальное оценивание математического ожидания ограниченной промежутком [-;] случайной величины с ПРВ f(x) как . Тогда ширина доверительного интервала для математического ожидания на основе оценок выборочного среднего и медианы - , в то время, как для ПИО-оценки (экстремальной порядковой статистики) , и это гарантирует большую эффективность ПИО-оценок для объемов выборки, начиная с некоторого достаточно большого N. Правда, на практике N может быть слишком велико (10³-10⁴ измерений).

В диссертационной работе доказано, что приведенные выше рассуждения имеют ограниченную область применимости: гарантируется только при f()>0, т.е. невырожденности помехи на границе диапазона. В противном случае , и ПИО-оценка может иметь меньшую эффективность по сравнению с альтернативными оценками медианы и выборочного среднего.

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

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

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

В рамках первой задачи получены оценки дисперсии медианной статистики в аналитической либо в рекуррентной форме: а) для равномерного закона, б) для симметричного распределения Лапласа, в) для распределения Симпсона и г) для обобщенно-степенной помехи с ПРВ при k>-1.

Вторая задача формализована следующим образом: выборка объемом N и медиана в окне (2l+1) измерение, позволяет работать с медианой как со случайной величиной, имеющей дисперсию , убывающую с расширением окна. Следует оценить, насколько выгодным является увеличение точности за счет расширения окна и повышения объема вычислений. Аналогом можно считать построение авторегрессионной модели по выборке N измерений с выбором оптимального порядка модели р [12]. Критериями является минимизация следующих функционалов качества:

финальная ошибка прогнозирования (ФОП) , (7а)

информационный критерий Акаике (ИКА) . (7б)

В работе обсуждена возможность использования тех же критериев для поиска оптимальной ширины медианного окна с формальной заменой pl.

В частности применение критериев (7а) - (7б) подтверждает: в случае распределений Лапласа и Симпсона процедура вычисления медианы по максимально большому числу отсчетов является выгодной, для обобщенного степенного закона с k=2 тенденция противоположная.

Известны оптимальные L-оценки для равномерного, лаплассовского и нормального распределения, но эта задача остается нерешенной в случае симпсоновской помехи. Заметим, что разностный сигнал от ЛТ1, аддитивно содержащего равномерно распределенную случайную помеху, представляет собой ЛТ0 с помехой Симпсона. Решена задача минимизации дисперсии линейной комбинации соответствующих 3-6 порядковых статистик и сформулированы общие закономерности формирования оптимальной L-оценки.

Цифровая фильтрация применительно к сигналам ЛТ0 и ЛТ1 - актуальная задача, т.к. аналитический пик при электрофорезе на микрочипах адекватно аппроксимируется кусочно-линейным сигналом, а параметры пика -амплитуда, временное положение и площадь пика, оцениваются с большими ошибками в условиях узости пиков или/и малых отношений сигнал/шум. Был использован тангенсный цифровой фильтр низких частот Баттерворта второго порядка, благодаря простоте реализации и достаточной эффективности подавления высокочастотных помех.

Рассчитано изменение формы ЛТ1 : передаточная характеристика фильтра (_1,2 - функции частоты среза, из условия W(1)=1, ), изображение исходного сигнала . Тогда изображение фильтрованного сигнала , а оригинал . Первое слагаемое есть x₁[k], а ₂ убывает с ростом относительной частоты среза f₀. Т.о., чем меньше частота среза, тем дольше продолжается колебательный переходный процесс. Т.о., возможна адекватная оценка параметров фильтрованного (деформированного) пика, благодаря подавлению шума, и последующее приведение параметров к исходному (дофильтрованному) пику.

Оценивание площади пика традиционно основано на численном интегрировании , где a и b - границы, F(x) - огибающая пика, f(x) - базовая линия. При малом отношении сигнал/шум, типичного для микроаналогов приборов электрофореза, сложность представляет как определение границ пика (а,b), так и адекватная оценка амплитуды сигнала, особенно, в области далекой от вершины, т.к. эти величины чувствительны к случайной составляющей. В работе исследовано упрощенное правило вычисления площади пика: , при h - высоте пика, _1/2 - ширина пика на половинной высоте. На примерах реальной электрофореграммы проиллюстрирована меньшая вариация упрощенной оценки площади.

Глава 4 посвящена анализу алгоритмов оценивания ЛТ0 и ЛТ1. В первой части представлено сравнение алгоритмов (рекуррентного и интервального) на информационной основе с констатацией очевидного вывода: в информациионном смысле рекуррентный алгоритм эффективнее интервального, особенно при большом динамическом диапазоне изменения сигнала.

Во второй части рассматривается алгоритм стохастической аппроксимации Роббинса-Монро (8) и его различные модификации:

(8),

где - оценка величины на n-ом шаге оценивания, y(c_n) - наблюдаемое значение величины y при , , есть искомый корень уравнения при М - математическом ожидании у. В дальнейшем процедура (8) получила развитие во многих работах: во-первых, разные подходы к выбору ; во-вторых, различный выбор функции ; в-третьих, модификации, позволяющие «отслеживать» точку на ЛТ1 (в частности, [3,4]).

Алгоритм Дупача свелся к следующему: имеются отсчеты измеряемой величины . Т.о. - детерминированный ЛТ1, - случайная помеха. Полагая - оценку на n-ном шаге оценивания, получим алгоритм:

(9).

Здесь - последовательность экстраполяций оценки. удовлетворяет условию несмещенности: при b=0 на любом шаге (т.е., ), при только асимптотически (при ). Параметр лежит в интервале (1/2;1). Интерпретировав (9) дискретной системой автоматического управления (САУ), получим, что корень передаточной характеристики линейной части будет равен , т.е. вне единичного круга на комплексной плоскости. Тем самым, линейная часть и все уравнение, безусловно, неустойчиво.

В диссертации предложена модифицированная экстраполяционная формула: . Т.к. , то несмещенность достигается на любом шаге. Передаточная характеристика линейной части соответствующей САУ имеет корень z=1 на границе области устойчивости, по крайней мере, кратности 2. В этом случае устойчивость негарантированна, но может достигаться за счет нелинейных слагаемых.

Тем самым, обоснована перспективность применения двухэтапного алгоритма оценивания параметра положения ЛТ1: переход к первому разностному сигналу (с возможной симметризацией помехи как свертки) и оценивание величины ЛТ0 с аддитивной помехой по сравнению с непосредственным отслеживанием ЛТ1.

Глава 5 посвящена анализу свойств алгоритма (8), главным образом, его модификации в форме Фабиана-Цыпкина. В первой части анализируется структура этого алгоритма и основные свойства оценки постоянного сигнала. Данная модификация имеет форму:

. (10)

При этом, , x_n+1=c*+x_n+1 -- измерение, x_n+1 -- помеха, с_n, c_n+1 -- оценки величины с* на n-ом и n+1-ом шагах оценивания, b -- параметр алгоритма, 2D -- величина зоны нечувствительности. В случае =0 , и получается сигнатурный алгоритм Цыпкина.

Анализ оценки алгоритма (10) предполагал исследование несмещенности, состоятельности и относительной эффективности в сравнении с оценкой МНК. Ранее [13] была доказана асимптотическая нормальность величины . Ее математическое ожидание равно нулю, а дисперсия . Исходя из необходимого условия экстремума , и минимальная дисперсия становится

.(11)

Отношение к дисперсии исходной помехи определяет эффективность оценивания. Это исследование было проведено ранее для некоторых помех (нормальная, равномерная, лаплассовская, смесь Тьюки). Новым явился анализ эффективности оценки при треугольной Симпсоновской помехе. Было впервые доказано, что (11) не зависит от величины . Причем, треугольный Симпсоновский закон единственный, обладающий таким свойством. Этот результат важен, т.к., в общем случае, можно оптимизировать подбором . Новый метод исследования сходимости оценки алгоритма Цыпкина предполагал преобразование (10) к САУ, что дало возможность использовать свойственные данной науке универсальные методы исследования.

Представленная САУ (см. рис. 2) характеризуется следующими свойствами: а) дискретна; б) имеет два контура обратной связи, в) имеет нелинейный элемент НЭ К₁ (неидеальное реле с зоной нечувствительности 2D, реализующее функцию y), г) содержит звено с переменными параметрами K₂[n], реализующее умножение на , д) линия задержки на 1 такт К₃ в цепи обратной связи. Соответственно x[n] -дискретная последовательность входных сигналов, с[n] - оценка c*, выполненная на n-ном шаге работы алгоритма (10).

Рис. 2. Структурная схема САУ, реализующей алгоритм (10) (пояснения в тексте).

Асимптотическая устойчивость положения равновесия САУ является необходимым условием несмещенности. Состоятельность оценки подразумевает ее асимптотическую сходимость к с* вне зависимости от входного воздействия, включая случайную помеху. Было осуществлено построение АЧХ всех элементов и линейной части САУ, и использован частотный критерий устойчивости Попова. В результате подтверждены следующие свойства оценки (10): а) при отсутствии зоны нечувствительности, т.е. при идеальном реле в качестве НЭ, САУ не имеет асимптотически устойчивого положения равновесия; б) в случае ограничения производной передаточной характеристики НЭ (в т.ч., ) асимптотически устойчивое положение равновесия достижимо; в) на устойчивость параметры помехи влияния не оказывают.

Другим оригинальным методом явилась модификация подхода Аоки [14] для анализа сходимости оценки (10). Для сходимости требовалось выполнение условия: при разность ПРВ оценок на (n+1)-м и n-м шагах 0 при любом . Используя алгебру событий, определение ПРВ и структуру (10) получим:

 (12)

Здесь ПРВ и функция распределения помехи , соответственно, f(x) и F(x).

Очевидно, что при требуемое условие выполнено. Допустив, симметричность помехи , можно на основе метода математической индукции по номеру шага n доказать симметричность любой ПРВ (12) .

Таким образом, модификация метода Аоки также позволяет утверждать, что оценка алгоритма (10) обладает состоятельностью и несмещенностью.

Глава 6 посвящена практической реализации алгоритма оценивания постоянного сигнала, что включает а) Подбор параметров (,), б) Выбор начального приближения оценки c₁, в) собственно оценивание, т.е. формирование последовательности c_n, г) формулирование критериев остановки оценивания (в т.ч., из-за разладки в последовательности измерений).

Решение первой задачи основывается на неравенстве . Т.е. необходима информация о р(х) - ПРВ помехи. Она решается введением начального участка длины M шагов, на котором рассчитывается размах помехи и средняя ПРВ. Если упорядочить выборку М измерений, можно получить оценку функции распределения помехи. Без упорядочивания, средняя ПРВ обратно пропорциональна размаху помехи - разности статистик x(1) и x(M), нахождение которых требует двух рекуррентных определений максимального и минимального значения.

Традиционный выбор оценки на первом шаге в форме обладает простотой; но отклонение c₁ от с* велико, если первое измерение - выброс, и алгоритм будет «настроен» на большую ошибку. Естественно выбрать с₁, основываясь на информации, уже полученной на начальном участке: а) при упорядочивании выборки, c₁ может быть комбинацией центральных статистик (медиана, оценка Гаствирта и т.п.), б) без упорядочивания можно выбрать с₁ как центр размаха . В случае выброса на первом измерении начальная ошибка оценивания уменьшится, по крайней мере вдвое.

Критерием остановки работы алгоритма (10) может явиться либо достижение требуемой точности оценивания, либо выявление разладки в последовательности измерений. Под разладкой понимается произошедшее на m-ом шаге скачкообразное изменение величины сигнала с* или его формы, например, из-за наложения линейного тренда.

При реализации алгоритма (10) вычисляется , и можно рассчитать -- сумму знаков поправок f_n в скользящем окне ширины L. Гистограмма распределения строится как рекурсивная процедура, что не требует значительных затрат. Показано, что при отсутствии разладки гистограмма симметрична относительно нуля. Количественно мерой асимметрии может быть аналог коэффициента асимметрии (нормированный начальный, а не центральный момент третьего порядка).

...