CC BY
60
УДК 681.3:519.9:519.68:612.001.57 А.М.Попов
РЕШЕНИЕ ОБРАТНОЙ ЗАДАЧИ ЭЛЕКТРОЭНЦЕФАЛОГРАФИИ С ПОМОЩЬЮ СТОХАСТИЧЕСКИХ МЕТОДОВ РАСПОЗНАВАНИЯ ОБРАЗОВ
(кафедра автоматизации научных исследований факультета ВМиК)
1. Введение. Работа посвящена развитию методов автоматизированного анализа электрической активности мозга. Электроэнцефалограммой [1] называется запись слабых (порядка 5-100 mV) электрических потенциалов, генерируемых мозгом. Сигнал представляет собой разность потенциалов между электродами, размещенными на поверхности головы. Еще Ханс Бергер [1] обратил внимание на то, что измеряемые пространственно-временные сигналы представляют некоторые образы, отражающие ощущения и состояние человека. Электроэнцефалография (ЭЭГ) — это метод, позволяющий судить о физиологической зрелости человека, его функциональном состоянии, наличии очаговых поражений, общемозговых расстройств и их характере. Одним из важных приложений ЭЭГ в медицине [2] является диагностика эпилепсии. Известно, что за генерацию эпилептического припадка отвечают локальные группы нервных клеток,так называемые очаги. На основе ЭЭГ-данных, а также поведения пациента во время приступа неврологи могут определить локализацию эпилептогенного очага.
Помимо медицинских приложений в последнее время одним из актуальных направлений в информатике является разработка так называемого интерфейса мозг-компьютер BCI (Brain-Computer Interface) [3, 4]. Появились работы, подтвердившие возможность достаточно достоверного распознавания образов на основе ЭЭГ, а именно выяснение вопроса, какую мыслительную деятельность из фиксированного набора возможных альтернатив развивает испытуемый [4].
Одной из важных задач автоматизированного анализа электрической активности мозга, возникающих из упомянутых приложений, является задача локализации источника (или зоны источников внутри мозга), создающего регистрируемый потенциал на поверхности головы. Задача локализации, рассматриваемая в настоящей работе, основана на дипольной модели источника [5-7]. Модель построена на предположениях о нейронных клетках как источниках (генераторах) электрического тока. Мозг моделируется как некий объемный проводник, обладающий неоднородной электропроводностью, с имплантированными в него диполями. Исходя из электропроводности тканей мозга можно показать,что время распространения сигнала очень мало, и математическая задача расчета потенциала сводится к решению краевой задачи для уравнения Пуассона [6, 7]. Таким образом, задача пространственной локализации источника решается на одном срезе по времени в ЭЭГ-сигнале. Пространственное решение в следующий момент времени не зависит от предыдущего момента. Однако при выполнении принципа непрерывности сигнала по времени можно предположить, что найденные пространственные локализации источника будут близки в близкие моменты времени.
Прямая задача ЭЭГ состоит в вычислении трехмерного электрического поля, создаваемого источниками-диполями [6, 7]. Полученное решение — потенциал электрического поля на поверхности головы — моделирует измеряемый в эксперименте потенциал электрического поля. Делая различные предположения о распределении и количестве диполей, можно добиваться совпадения экспериментального потенциала с модельным. Прямая задача должна учитывать структуры тканей мозга с различными проводимостями и другие геометрические параметры для правильного вычисления индуцированного потенциала.
Обратная задача ЭЭГ состоит в определении положений источников (правой части уравнения Пуассона) по измеренным потенциалам на поверхности головы [8-10]. Для решения обратной задачи требуется решение большого числа прямых задач, эффективному решению которых уделяется много внимания. Известно, что существует бесконечное число конфигураций внутренних источников тока, которые могут создавать любое заданное на поверхности головы распределение потенциала [11]. Это означает, что обратная задача ЭЭГ в общей постановке не имеет единственного решения. Чтобы сделать решение задачи единственным, приходиться вводить некоторую модель источников тока. Одной из таких моделей является модель дипольных источников. Во многих случаях используется один (эквивалентный) диполь. Неизвестных в задаче в этом случае шесть: три пространственные координаты
расположения диполя и три компоненты вектора тока, характеризующего направление и величину момента диполя.
Существует два подхода к локализации источников. Первый подход основан на фиксированном расположении диполей, а задача сводится к линейной, т.е. к отысканию их силы тока (задача компьютерной томографии) [12]. В этом случае задача решается на основе методов регуляризации, предложенных в [13, 14]. Второй подход основан на поиске как положения диполя, так и его момента ("подгонка диполя") [10, 15]. В этом случае задача нелинейна. В настоящей работе именно в этой постановке рассматривается задача локализации источника.
Учет временной составляющей сигнала может уточнить локализацию источника, если учесть статистику появления источника в различных областях. Модели анализа пространственно-временного сигнала ЭЭГ можно найти в [16, 17]. В модели "изучения отклонений" (deviation scan) предполагается подгонка одиночного диполя, повторенная для большого числа его расположений. Модель многократной классификации сигнала MUSIC (Multiple Signal Classification) [6] определяет алгоритм поиска источников, являющихся независимыми друг от друга во времени. В модели "распределенных источников" предполагается, что источники распределены по некоторой заранее известной поверхности. Часть поверхности коры головного мозга, на которой расположены источники, определяется априори. В работе [17] рассматривалась модель подгонки нескольких дипольных источников. Основная идея работы состоит в использовании Генетического алгоритма для ЭЭГ-сигнала при пространственной реконструкции источников.
В настоящей работе анализ сигналов ЭЭГ проводится на основе метода Random Forest [18]. Алгоритм Random Forest является одним из наиболее многосторонних алгоритмов классификации данных, известных в области Data Mining. Дипольная функция источника позволяет моделировать электрический потенциал, измеренный на поверхности головы, через параметры дипольного источника — его координаты (пространственную локализацию) и вектор момента. Область наиболее активных источников восстанавливается по экспериментальным данным в выбранных срезах по времени. В работе предлагается алгоритм вычисления границ области локализации, основанный на аппарате деревьев решений [19]. Для каждого момента времени несколько диполей случайным образом располагаются внутри области, и строится свое дерево решений, классифицирующее набор параметров дипольных источников на два класса по уровню допустимой ошибки приближения ЭЭГ-потенциалов. Каждое дерево из временного ансамбля оценивает область локализации дипольных источников, вычисляя порог классификации. К первому классу относятся диполи, дающие ошибку в восстановлении потенциала меньше задаваемого порога. Классификация проходит без усечения ветвей на случайно выбранном для каждого момента наборе признаков — параметров диполя. Построенный таким образом ансамбль деревьев классификации представляет собой случайный лес [18]. Решение о зоне локализации принимается, если оно подтверждено большинством деревьев из леса.
2. Дипольная модель источников. В дипольной модели мозг рассматривается как объемный трехмерный проводник. Обозначим VQi его объем и Sioi — поверхность, ограничивающую объем. Источниками электрической активности являются электролитические токи внутри нервных клеток коры головного мозга [15]. Закон Ома в рассматриваемом проводнике запишется в следующем виде:
j=jm + o--E, (1)
где jin — плотность стороннего (ионного, электролитического по происхождению) тока, создаваемого нервной клеткой. Источники тока помещены в проводящую среду с неоднородной электронной проводимостью (7(г). Ионные токи в объеме мозга порождают электрическое E(r, i) и магнитное B(r, i) поля, которые описываются уравнениями Максвелла.
Показано, что в данной задаче временная и пространственная составляющие уравнений могут быть разделены и пространственная часть в каждый момент времени удовлетворяет стационарному уравнению
(V-j) = 0. (2)
Применяя оператор дивергенции к обеим частям (1), используя (2) и вводя потенциал электрического поля U(г, t):
Е = -VU, [V х Е] = 0, (3)
получаем уравнение Пуассона:
(V • <70VC/((r, t))) = —(V ■ iin). (4)
Предполагая, что все источники тока заключены в объем VQi, а в остальном пространстве Vrest источников нет, можно записать частное решение уравнения Пуассона в виде объемного потенциала в
Области С ПОСТОЯННОЙ ПРОВОДИМОСТЬЮ (То-
Задача состоит в вычислении потенциала U как решения уравнения (4). Для этого потенциал U представляется в виде суммы:
U(r1t) = W(r1t) + V(r1t)1 (5)
где W — потенциал, создаваемый диполями:
= (6)
0 RMPi J
Здесь Ni — число диполей и V(r,t) — неизвестный потенциал индуцированного поля, который создается из-за наличия границ и неоднородной электропроводности. Если известны положения диполей и их моменты (прямая задача), тогда потенциал W(rM,t), создаваемый диполями, известен и требуется найти индуцированный потенциал V(r,t). В областях однородности проводимости для функции V(r,t) имеем задачу Неймана для однородного уравнения Лапласа
AV(r,i) = 0 (7)
с неоднородным граничным условием на скальпе (граница с непроводящей средой):
дУ _ dW дп дп
и условиями сшивки на каждой поверхности Sk разрыва проводимости:
(8)
Vk = Vk+1
(9)
Sk
Заметим, что непрерывно вместе со своими производными во всем пространстве за исключением точек расположения источников.
Нормальные производные потенциала V терпят разрыв:
dvk dvk+1 dw
= °к+1 • —+ К+1 - ы
дп дп дп
(10)
Sk
Здесь индекс к обозначает номер поверхности, ограничивающей слой проводимости. При учете проводимости области вне головы (проводимость геля электродов, учета других проводящих слоев вне скальпа) задача ставится в неограниченной области.
Метод решения прямой задачи основан на использовании общего решения уравнения Лапласа в виде сферических гармоник [7]. Коэффициенты общего решения находятся методом наименьших квадратов для удовлетворения граничного условия Неймана на произвольной (не сферической) границе и условий сшивки на границах смены электропроводности. В общем случае будем считать, что к поверхностей, к = 1, 2, 3,..., К, заданных уравнениями = ($, ф) в сферических координатах, разделяют области с различными проводимостями. Здесь к = 1 обозначает самую внутреннюю поверхность, содержащую начало координат, и к = К — поверхность скальпа (внешняя). Пусть Q = {($¿, (fij) I i = 1,..., Ntf] j = 1,..., iV^} — сетка по ■& и ср. Обозначим = ($¿, tpj).
Ищем решение в к-й области в виде разложения по сферическим гармоникам:
Г "
Vk(r, i}, р) = | г" (Акпт cos т<р + Вкпт sin (cos #) +
п=0 ш=0
П ч
+ г-(п+1) ^^ {Скт COS т<р + Dkm sin ТП(р)р(т^ (cos ■&) У. (11) n J
m = 0
Для любых коэффициентов А, В, С ж И функции Ук удовлетворяют уравнению Лапласа (7). Коэффициенты Акт1 Вкт1 Скт, Окт, к = 1,..., К, находятся из условий сшивки (9), (10) на произвольной
(к)
поверхности г\- .
Для к = 1 коэффициенты С„т = 0]1т = 0 из требования ограниченности в нуле, и для к = К = = 0 из ограниченности на бесконечности. Остальные коэффициенты находятся из условий сшивки (9), (10) в смысле наименьших квадратов:
[Vk(rij(k),-&i, <~Pj) - Vk+1(rij(k),{>i, cpj)]2 min,
J
£2=5]
gyk dVk+1
yij{k),"&i,(pj) ~ Vk
dW
dn
dn (rij(k),#i, 4>j) ~ (fij)
mm,
(12) (13)
где Vk = z%L
Необходимые условия минимума функционалов приводят к системе линейных алгебраических уравнений относительно которая решается с использованием SVD (Singular
Value Decomposition) алгоритма [7].
3. Задача локализации нейронного источника. Математическая постановка обратной задачи сводится к нахождению правой части неоднородного трехмерного эллиптического уравнения (4) с неоднородными коэффициентами по решению С/, заданному в точках поверхности 5теав, и граничному условию Неймана на этой поверхности.
Для численного решения обратной задачи используются приближенные граничные условия в виде условий минимума функционалов ошибки:
ef>, гр) = \\Uexp - W(u, гр, г) - Vk(t) d[W{u, rpir) + VK(r)]
дп
mm,
mm .
(14)
(15)
Задача состоит в определении величин и координат диполей V, гр для наилучшего приближения потенциала итос1е1 = ^(и, гр, г) +Ук(г) на поверхности скальпа к экспериментально измеренному распределению потенциала иехр путем минимизации функционала е = е^ + е^. Для итеративного решения обратной задачи необходимо решать прямую задачу на каждой итерации, т.е. по дипольным источникам с известными расположением и мощностью находить трехмерное распределение потенциала.
4. Алгоритм классификации источников с помощью деревьев решений. Анализируемая электроэнцефалограмма представляет собой зависимость измеренного потенциала от времени на 14-ти сенсорах поверхности головы. Для локализации диполя мы рассматриваем лишь небольшое временное окно от £ = 1,2 с до £ = 1,3 с. В этом окне 10 временных точек. Обучение на тренировочном наборе данных соответствует известному алгоритму С4.5 [19].
Пусть нам задано некоторое обучающее множество Т, содержащее образы, каждый из которых характеризуется М признаками (атрибутами). Пусть через С^Сг, • • •, С& обозначены метки классов. Обозначим вектор признаков через
ХР = {Х1,Х2, • • .,Хг, . ..,Хм\Ск}Р,
где р — номер примера, к — номер класса, х? — г-й признак.
Случайным образом в момент времени Ь1 располагаются Мр (полное число примеров) диполей внутри головы (полной допустимой области) и выбираются их моменты. Здесь (х^х^х3,...) параметры диполей:
х1 = Ы} = {г*,^,^,^,^,^},
2 г2т г2о2 2 2 2 2т Х = {*к} = {Гр,др,<рр,1У
х3 = {х3,} = {г3,*д3 ч?3 и3,и3д, V3},
Каждый из наборов параметров (х1, х2, х3,...) принадлежит соответственно первому (х1), второму (х2), третьему (х3) и т.д. диполям. Для каждого р-го диполя вычисляется ошибка
£Р = ~ - I sin eihtihyj. (16)
i 3
Задаем уровень допустимой ошибки eth и разделяем параметры диполей на два класса: ниже порога, когда ер < Eth (Ci — первый класс), и выше порога (С2 — второй класс). Таким образом, получен тренировочный набор данных. Основой классификации является уровень допустимой ошибки по потенциалу. Этот уровень определяет точность классификации. На рис. 1 ,a,b,c,d представлены линии уровня ошибки по потенциалу е = RRE на плоскости двух параметров — координат диполя (гр — вр). Линии уровня представлены для четырех моментов времени: t = 0,1,2,3. Здесь использовано международное обозначение ошибки по потенциалу RRE (Residual Ralative Error). Результаты, представленные на рисунке, позволяют судить о связи ошибки по потенциалу с ошибкой в локализации диполя. Конечно, это суждение здесь очень приблизительно, так как представлена только одна проекция в шестимерном пространстве параметров. Хорошо видна сложность задачи локализации — это наличие очень плоских оврагов в функционале ошибки. Другой принципиальной чертой задачи распознавания является ее неустойчивость по отношению к тестовым данным. Чем меньше ошибка по потенциалу, тем больше вероятность того, что области локализации по различным деревьям не перекрываются. В этом случае вероятность CTDVKTVDHoro риска повышается.
в г
Рис. 1. Линии уровня ошибки при аппроксимации потенциала RRE = const на плоскости координат диполя (Гр-вр) для четырех последовательных моментов времени: t = 0 (a), t = 1 (б), t = 2 (в), t = 3 (г)
Узлом дерева является узел проверки, а листом — узел решения. Условие в узле будем называть тестом. Множество Т содержит примеры, относящиеся к разным классам. Цель построения дерева — разбить множество Т на некоторые подмножества Т{, которые будут состоять в основном из примеров,
относящихся к одному классу. Нас интересует разбиение внутренней области на подобласти, в которых примеры относятся к первому классу, т.е. такому, в котором ошибка по потенциалу наименьшая.
Вначале выбирается наиболее значимый атрибут для разбиения. Выбранный атрибут должен разбить множество так, чтобы получаемые подмножества состояли в основном из объектов, принадлежащих к одному классу. Количество образов-примесей в каждом из этих множеств будет наименьшим. Мы используем теоретико-информационный критерий для выбора наиболее значимого атрибута. Для различных классов и различных примеров каждый признак хг упорядочивается в порядке возрастания величины. Вычисляются средние значения соседних величин признака:
X* = (Х, + Х, + 1)/2.
Эти значения х* являются границами разбиения на подмножества Т{ по признаку ж^,г = 1,...,М—1. Пусть N¡1,1 есть число примеров из Т^, относящихся к одному классу С к- Тогда вероятность того, что случайно выбранный пример из множества Т{ будет принадлежать к классу Ск, равна = А^/А^, где Ni есть полное число образов в множестве Информация, содержащаяся в выборе равна = — 1с^2 Рк,г, а средняя информация, приходящаяся на такой выбор, есть энтропия:
" А-
Н(хг, Г,) = - ^{Лм • 1оё2 Рк,г} • (17)
к= 1
Обозначим через энтропию полного множества Т при анализе по признаку Х{. Тогда крите-
рием выбора признака будет минимальное значение
С(А) = Н{хг,Т) -Н{хг,Тг) (18)
по всем признакам. Итак, множества ... ,Тп получены при разбиении исходного множества Т
по признаку Х{. Выбирается атрибут, дающий максимальное значение по критерию (17).
Пусть выбран признак хг для всех примеров. Теперь разбиваем множество Т на два подмножества путем вычисления порога по выбранному признаку. Для этого вычисляется пороговое значение признака ж,^/, в наборе ж* из условия максимизации энтропии (17), т.е. тахЯ(х*, Т^. Исходное множество
' X*
г
примеров разбивается на два подмножества Т^ и Т^. Множество Т^ соответствует критерию
Хг < Хг,ш, (19)
а множество Т^ определяется неравенством
Хг > Хг,м- (20)
Далее рассматриваем множество Т^. Делаем проверку: если в Т^ находятся примеры одного класса, то получено решение и алгоритм останавливается; если в Т^ присутствуют примеры других классов, то продолжаем разбиение каждого из полученных множеств Т^ и Т^. Затем анализируем дальше Т^, в котором выбираем наиболее информативный признак для следующего разбиения, используя указанную выше процедуру. Причем признак Xi исключается из соревнования. Пусть выбран признак х^ для следующего разбиения. Из условия максимизации энтропии вычисляется порог х^^ь, по признаку х^ и множество Т^ оказывается разбитым на два: Тьь и ТЬК по признаку х^. Если полного разделения образов не получилось, то разбивается множество Т^ на два множества: ТКЬ и ТКК. Таким образом, используя два признака, мы разбили исходное множество на 4 подмножества: Тьь, ТЬК, ТКЬ и ТКК. Далее процедура повторяется. В данной работе используется определенный критерий остановки алгоритма (ограничение глубины дерева). Когда проведено разделение множества примеров по всем признакам, алгоритм останавливается. После остановки алгоритма принимается решение об отнесении образов к классу и вычисляется вероятность принадлежности образа к классу.
Классификация новых данных происходит с помощью дерева, и вычисляется ошибка классификации. Пусть Ак — число примеров из ТЬК, относящихся к одному классу Ск- Тогда вероятность того, что случайно выбранный пример из множества ТЬК будет принадлежать к классу С к, равна Рк = Nk/NzJR1 где А^д есть полное число образов в множестве ТЬК. Если Рк > где — заранее выбранное пороговое значение вероятности классификации. Например, если число признаков равно двум, то дерево содержит три узла (пороговые значения для разбиения) и ссылку на потомков. Если
остановки не произошло по ходу разбиения, то дерево содержит четыре листа, соответствующие решениям, принятым в четырех подмножествах множества образов. Эти пороги xth, указатели и листы (Pth) должны содержаться в памяти для классификации образов.
5. Локализация положения диполей методом Random Forest. В случайном лесе выращивают много деревьев классификации. Для того чтобы классифицировать новый образ из входного вектора, входной вектор подается на каждое дерево из леса. Каждое дерево дает классификацию, и мы можем говорить о том, что оно голосует за определенный класс. Лес выбирает классификацию, имеющую наибольшее количество голосов по всем деревьям из леса.
Каждое дерево строится следующим образом.
1. Если число примеров в тренировочном наборе есть Np, то выбираем Np примеров случайным образом. Эта выборка будет тренировочным набором данных для выращивания дерева.
2. Если число признаков образов есть М (в нашем случае шесть), то задается число то <С М такое, что на каждом узле то признаков выбираются случайным образом из М и их число поддерживается постоянным в процессе выращивания леса.
3. Каждое дерево растет до наибольшего возможного размера.
Покажем работу алгоритма на тренировочной базе данных. Предлагаемый алгоритм локализации, т.е. нахождения возможных координат диполя, следующий.
1. Случайным образом в момент времени t1 располагаются Np диполей внутри головы (полной допустимой области) и выбираются их моменты.
2. Для каждого р-го диполя вычисляется ошибка
£Р = ~ - I #i,<Pj)]2 sin tiihtihvj. (21)
i 3
3. Задаем уровень допустимой ошибки eth и разделяем базу данных (параметры диполей) на два класса: ниже порога (первый класс) ек < eth и выше порога (второй класс). Таким образом получен тренировочный набор данных.
4. Для классификации этих данных применяем дерево решений, но используя лишь три (то = 3) случайно выбранных признака. Так как набор признаков мал, то это дерево строится без усечения. Выходом построения первого дерева будут три параметра (г^, $fh, <~p\h) — пороги, определяющие трехмерную область локализации диполя первого класса.
5. Для второго момента времени процедура повторяется и строится второе дерево из леса. Выходом построения второго дерева будут свои три параметра — пороги (rjh, i)2th, определяющие трехмерную область локализации диполя первого класса.
6. Процесс построения леса повторяется и выбирается такая классификация (локализация диполя), за которую проголосовало большинство деревьев.
На рис. 2 показана последовательность зон локализации источника, соответствующая четырем моментам времени. Видно, что в процессе обучения на тренировочном наборе данных область перекрытия зон, найденных деревьями, становится уже и искомое положение диполя находится в этой зоне. Мы можем говорить о сходимости к некоторому среднему положению диполя.
Качество использованного алгоритма классификации зависит от задаваемой точности восстановления потенциала. Поэтому необходимо определить зависимость точности определения локализации диполя от задаваемой погрешности ек по потенциалу. Ошибка локализации диполя есть сумма ошибки на тренировочных данных и вероятной ошибки на тестовых данных. Чем меньше ошибка по потенциалу, тем уже области локализации в каждый момент времени, тем более вероятно, что зоны перекрытия у большинства деревьев не будет и повышается вероятность ошибки на тестовых данных. Для уравновешивания этих ошибок вводится параметр S регуляризации по областям локализации:
£ — Strain + 5 (22)
Чем больше S, тем шире зона локализации и меньше точность. Кривые, приведенные на рис. 3, показывают, что при небольшой погрешности по потенциалу, меньше 5%, погрешность локализации диполя становится того же порядка. При большой погрешности по потенциалу сходимость теряется и окончательная зона перекрытия становится слишком большой, чтобы можно было говорить о локализации. Дальнейшее уменьшение задаваемой погрешности было бы некорректно.
Рис. 2. Области локализации диполя (на плоскости гр—вр), определенные деревьями решений для четырех последовательных моментов времени: а, Ь, с, г1. Выделенное пересечение прямоугольников соответствует области
локализации источников
Рис. 3. Зависимость ошибки локализации диполя БЬЕ от ошибки аппроксимации потенциала на измерительной поверхности РШЕ: кривая 1 соответствует параметру регуляризации 5 = 0,01, кривая 2 — параметру
регуляризации 5 = 0,04
В [18] показано,что величина ошибки определяется двумя причинами: во-первых, корреляцией между любыми двумя деревьями в лесу; во-вторых, силой каждого индивидуального дерева в лесу. Возрастание корреляции приводит к возрастанию ошибки леса. Дерево с малой величиной ошибки является сильным (хорошим классификатором). Возрастание силы индивидуальных деревьев уменьшает ошибку леса.
Уменьшение числа т уменьшает как корреляцию, так и силу дерева. Увеличение т приводит к увеличению обоих этих факторов. Где-то между ними находится оптимальная область значений т. Обычно это достаточно широкая область. Используя известную оценку ошибки, значение т в области может быть быстро найдено. Это единственный подгоночный параметр, по отношению к которому случайные леса чувствительны.
Мы провели численное исследование с целью определения оптимального числа случайно выбираемых признаков из полного набора 6-ти признаков. В расчетах число т было выбрано равным трем. При этом лучшая достигнутая точность локализации составила 2,1%.
6. Заключение. В настоящей работе предложен алгоритм анализа сигналов ЭЭГ на основе метода Random Forest [18]. В рамках дипольной модели источников электрический потенциал, измеренный по нескольким каналам на поверхности головы, выражен через параметры дипольного источника — его координаты (пространственную локализацию) и вектор момента, которые составляют пространство признаков задачи. Зона локализации источника является одним из выходов деревьев классифи-
кации в каждом срезе по времени. Ансамбль деревьев решений ищет перекрытие зон локализации, полученных каждым деревом в отдельные моменты времени. Показано, что использование выборки признаков, равной трем в каждом дереве, наиболее эффективно и устойчиво решает проблему локализации, при этом величина структурного риска классификации на тестовых данных минимальна. Относительная погрешность определения зоны локализации — того же порядка, как и априори задаваемая погрешность в восстановлении потенциала. Предложенный алгоритм требует одного решения прямой задачи ЭЭГ в каждый момент времени. Вычисления могут быть ускорены за счет параллелизма самого алгоритма Random Forest. Приведенные экспериментальные расчеты с модельной базой данных, допускающей аналитическое рассмотрение, показывают эффективность предложенного подхода. Алгоритм хранит данные классификации в наиболее сжатом виде,что позволяет классифицировать данные ЭЭГ в реальном времени.
Автор благодарит чл.-корр. РАН JI.H. Королева за обсуждение результатов работы и сделанные замечания и Е.А. Попову за создание программного комплекса Random Forest и проведение расчетов.
СПИСОК ЛИТЕРАТУРЫ
1. Tripp J., Physical concepts and mathematical models // Biomagnetism: An Interdisciplinary Approach / Eds.: S.J. Williamson, G.L. Romani, L. Kaufman, I. Modena. N.Y.: Plenum, 1983. P. 101-139.
2. ГнездицкийВ.В. Обратная задача ЭЭГ и клиническая электроэнцефалография. Таганрог: Изд-во ТРТУ, 2000.
3. Vi dal J.I. Real-time detection of brain events in EEG // Proc. IEEE. May 1977. 65. P. 633-664.
4. Wolpaw J. R., В i r b a u i n e r N.,McFariand D. J., Pf u r t sc h e 11 e r G.,Vaughan T.M.Brain-computer interfaces for communication and control // Clin. Neurophysiol. 2002. 113. P. 767-791.
5. Plonsey R., Heppner D.B. Considerations of quasistationarity in electrophysiological systems // Bulletine of Mathematical Biophisics. 1967. 29(4). P. 657-664.
6. Mosher J. C., Leahy R. M., Lewis P. S. EEG and MEG: Forward solutions for inverse methods// IEEE Trans. Biomed. Eng. March 1999. 46. N 3. P. 245-259.
7. Хоффманн К., Попов A.M., Певцов C.E., Федулова И.А. Численное решение прямых математических задач электроэнцефалографии // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2004. № 3. С. 19-26.
8. Захаров Е. В., К о п т е л о в Ю. М. О некоторых математических проблемах в решении обратной задачи электроэнцефалографии // ДАН СССР. 1987. 292. С. 578-581.
9. SarvasJ. Basic mathematical and electromagnetic concepts of the bio-magnetic inverse problem // Phys. Med. Biol. 1987. 32. P. 11-22.
10. Хоффманн К., Попов A.M., Федулова И. А., Певцов С.Е. Численное решение обратных математических задач электроэнцэфалографии // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2004. № 4. С. 16-27.
11. von Helmholtz H.L.F. Uber einige Gesetze der Verteilung elektrischer Strome in körperlichen Leitern mit Anwendung auf die tierischelektrischen Versuche // Ann. Physik und Chemie. 1853. 89. S. 211-233, 354-377.
12. Pascual-Marqui R. D. Review of methods for solving the eeg inverse problem // International journal of bioelectromagnetism. 1999. 1(1). P. 75-86.
13. Тихонов A.H., Арсенин В.Я. Метоы решения некорректных задач. М.: Наука, 1986.
14. Наттерер Ф. Математические аспекты компьютерной томографии. М.: Мир, 1990.
15. Wolters С.Н. et al. Comparing regularized and non-regularized nonlinear diole fit methods: A study in a simulated sulcus structure // Brain Top. 1999. 12(1). P. 3-18.
16. Schmitt U., Louis A.K. Efficient algorithms for the regularization of dynamic inverse problems. Part I: Theory // Inverse Problems. 2002. 18. N 3. P. 645-658.
17. Хоффманн К., Попов A.M., Федулова И.А., Певцов С.Е. Моделирование пространственно-временной электрической активности нейронных источников // Прикладная математика и информатика / Под ред. чл.-корр. РАН Д.П. Костомарова и ак. РАЕН В.И. Дмитриева. М.: Макс-Пресс, 2004. № 17. С. 55-71.
18. В г ei man L. Random Forests // Machine Learning. 2001. 45. P. 5-32.
19. Ross Quinlan J. C4.5: Programs for machine learning. Kaufmann Publishers, 1993.
Поступила в редакцию 10.02.06
УДК 517.968
В.И.Дмитриев, Е.В.Захаров
ГРАНИЧНЫЕ ИНТЕГРАЛЬНЫЕ УРАВНЕНИЯ ПЕРВОГО РОДА В ЗАДАЧАХ ЭЛЕКТРОДИНАМИКИ НЕОДНОРОДНЫХ СЛОИСТЫХ СРЕД
(кафедра математической физики факультета ВМиК)
Метод граничных интегральных уравнений Фредгольма второго рода является основным математическим аппаратом теоретического исследования различных задач для уравнений в частных производных. Классическая теория потенциала и разнообразные ее обобщения позволяют переходить от дифференциальных уравнений к интегральным, а затем изучать проблему разрешимости.
Построение схем и разработка вычислительных технологий для численного решения интегральных уравнений были начаты с конца 50-х гг. и естественно были связаны с решением интегральных уравнений Фредгольма второго рода. Интегральные уравнения первого рода либо вообще не рассматривались, либо служили основой построения различных асимптотик и точных аналитических решений.
В связи с созданием теории некорректно поставленных задач (к которым относятся и многие интегральные уравнения первого рода) и разработкой методов их численного решения возникла необходимость создания устойчивых вычислительных технологий решения интегральных уравнений первого рода, непосредственно связанных с краевыми задачами электродинамики. В середине 60-х гг. был открыт и освоен метод численного решения интегральных уравнений Фредгольма первого рода с логарифмической особенностью при совпадении аргументов [1-3]. Был преодолен психологический барьер невозможности применения таких уравнений и началось интенсивное использование интегральных уравнений первого рода для разработки вычислительных технологий исследования электромагнитных полей в неоднородных средах, в том числе и слоистых.
В данной публикации рассмотрены некоторые применения метода граничных интегральных уравнений в задачах электродинамики неоднородных слоистых сред. При этом без ограничения общности рассмотрение проведено для плоских задач электродинамики неоднородных слоистых сред.
1. Дифракция электромагнитного поля на идеально проводящей поверхности, расположенной в слоистой среде. Пусть имеется плоскопараллельная слоистая среда, характеризуемая кусочно-постоянными комплексно-значными функциями е(г) и В данной среде в одном из сло-
ев расположена идеально проводящая бесконечно тонкая цилиндрическая поверхность 5, бесконечно протяженная вдоль оси Ох. Обозначим Г контур (замкнутый или разомкнутый) в плоскости х = 0 и будем считать контур Г достаточно гладким.
Рассмотрим гармоническое во времени (е-гс^) электромагнитное поле {Е, Н}, возбуждаемое в среде произвольной системой источников, плотность распределения которых не зависит от х. Нормальное поле в слоистои среде обозначим {Е°,Н0}.
Могут быть поставлены следующие краевые задачи на плоскости уОг.
1. Е-поляризованное поле Е{ЕХ, 0,0}, Н{0,Уу,0}
