Анализ соотношений главной последовательности САККАД с позиции оптимального быстродействия

Саккадическое движение глаза относится к специфическим быстрым движениям человека и играет ключевую роль в повседневной активности. Саккады анализируются во многих областях исследований, например, в психологии, медицине, спорте [1–3].

В последние годы достигнут значительный прогресс в понимании механизмов, лежащих в основе саккадических движений, что стало возможным благодаря развитию технологий регистрации движений глаз и широкому применению математического моделирования. Особое внимание при исследовании саккад уделяется саккадической главной последовательности (главной последовательности) – зависимости между амплитудой саккады, ее максимальной скоростью и длительностью [4; 5]. В литературе предложен ряд параметрических моделей и методов измерений, направленных на точное описание динамики движений глаз [6; 7].

Проблемой имеющихся описаний главной последо-

0009-0005-4816-3532

0000-0001-9720-8163

Эта статья доступна в соответствии с условиями лицензии Creative Commons Attribution-NonCommercial 4.0 International

License (CC BY-NC 4.0)

This work is licensed under a Creative Commons Attribution-NonCommercial 4.0 International License (CC BY-NC 4.0)

вательности (ГП) для саккад и самих саккад следует считать вариативность применяемых математических описаний и критериев, зачастую без всякого обоснования выбора модели. В подобной ситуации представляется рациональным обратиться к эталонным саккадам. Для формирования идеальной, модельной саккады используются модели другого рода. Зачастую они основываются на физиологических и механических представлениях о работе системы [8–11].

Одними из наиболее обоснованных представляются модели, использующие механическую модель глаза и рассматривающие саккаду как оптимальное по времени движение. На практике системы управления движением глаз находят всё более широкое применение в современных технологиях человеко-компьютерного взаимодействия, включая интерфейсы и системы виртуальной реальности [12–16]. Разработка таких систем требует использования математических моделей как для отдельных саккадических движений, так и для параметров глазодвигательной активности, к которым относится ГП.

Движения глаз обеспечивают возможность стабилизировать и переводить направления зрительной оси глазного яблока в пространстве. Саккады, произведенные в ответ на предъявленный стимул, – это быстрые и точные повороты глаз, обеспечивающие перевод взора от одного объекта к другому. Эти движения обеспечиваются глазодвигательным аппаратом и структурами нервной системы и тесно связаны с визуальным восприятием, вниманием и характеристиками окружающей среды. Считается, что саккада – рефлекторное движение, при котором нервная система решает задачу наискорейшего перевода взора (поворота глаза). Таким образом, представляется целесообразным моделировать ее задачей быстродействия [8; 10].

Следуя [10], опишем глазное яблоко как твердое тело, будем учитывать влияние вязкого сопротивления и упругих сил, действующих со стороны глазодвигательных мышц, зрительного нерва и окружающих тканей. Будем рассматривать одномерное движение в аксиальной плоскости. Введем обозначения для параметров глазного яблока: R – радиус, I – момент инерции, φ – угол поворота; ф - соответственно угловая скорость; и параметры окружающих тканей и глазодвигательных мышц, реализующих движения глаза: ν – коэффициент вязкости, k - коэффициент жесткости, X 1 , X 2 - управляющие параметры мышц агониста и антагониста [17; 18], u ₁, u ₂ – скорости сокращения мышц (величины наиболее простым образом описывающие изменения

X 1 , X 2 - управляющие сигналы). Система, описывающая динамику глаза, как и в [10] имеет вид:

JIф + 2vR2ф + 2kR2ф = kR (Xj - X2),

[X. = u. , i = 1,2.                                        ⁽¹⁾

Вводя обозначения для переменных системы x1 = ф, x2 = ф, x3 = k (X1 - X 2) /1, управление

u = u 1 - u 2 , и коэффициенты A 1 = 2 vR ² / 1 , A 2 = 2 kR ² / 1 , из (1) получим систему в форме Коши:

x 1 = ^x 2

‘ x 2 =^- A 2 x i ^- Ax 2 + x 3 .                 (2)

x ₃ = u

Граничные условия . Начальные условия без ограничения      общности      выберем      нулевыми:

x_x ( 0 ) = x ₂ ( 0 ) = x ₃ ( 0 ) = 0, в конечный момент времени T :

x j ( T ) = A - амплитуда саккады, x ₂ ( T ) = 0, x ₃ ( T ) = Л_Т .

x3 имеет размерность удельного момента силы. Очевидно, что как сила мышцы, так и скорость ее изменения – величины ограниченные. Формализуем это как ограничение на управление, некоторой величиной u . Тогда построим u минимизирующий функционал:

J = T ^ min, | u | <  u * .

u

Воспользуемся, как и в [10], обезразмериванием системы и перейдем к переменным x. * y. = x., i = 1,3, u = u*v с размерными коэффициентами x.*. Введем характерные величины времени и амплитуды T* = 0,1 c, x1* = n / 9, остальные величины и коэффици енты выражаются через них:  v* = x3* / T*3, A1T* = a1,

A 2 T * ² = a ₂, a 1 = 2,71, a ₂ = 1,43.

Во введенных обозначениях задача быстродействия примет вид:

J y ¹ = y ²

⁴ y 2 =^- a 2 У 1 ^{- a} 1 У 2 ⁺ У 3 ,

. y 3 = v

У 1 ( ⁰ ) = У 2 ( ⁰ ) = У 3 ( ⁰ ) = ⁰,                  ⁽⁴⁾

У 1 ( T ) = A t , У 2 ( T ) = 0, У 3 ( T ) = Л T ,

J = T ^ min, | v | <  v * .

v

На основе принципа максимума Понтрягина [19] и теоремы Фельдбаума [20; 21] управление системой (4) предполагается кусочно-непрерывной функцией следующего вида:

v*,  0 < t < t1, v = ^ -v*,  t1 < t < 12,

Kv*,  12 < t < T, где t1 , t2 – времена переключений управления.

Собственные числа системы (3):   X 1 « 0,718,

а

б

Рис. 1. Изменения углов поворота ( a ) и скорости глаза ( б ) для саккад амплитудой 6,8°, 10°, 15°, 20°, 30°, 40°

Таблица 1

Длительности интервалов постоянного управления для выбранных амплитуд саккад, % от продолжительности соответствующей саккады

A	6,8°	10°	15°	20°	30°	40°
Этап 1
Случай 1	29,37	30,71	31,50	32,26	33,53	34,29
Случай 2	28,92	30,52	31,48	32,25	33,45	34,18
Этап 2
Случай 1	49,28	50,00	50,00	50,00	50,00	49,53
Случай 2	48,60	49,8	50,00	50,00	50,00	49,32
Этап 3
Случай 1	21,35	19,29	18,50	17,74	16,48	16,18
Случай 2	22,49	19,68	18,53	17,76	16,55	16,50

Таблица 2

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

T	6,8°	10°	15°	20°	30°	40°
Случай 1	49,97	55,11	63,36	69,60	79,43	88,02
Случай 2	50,54	55,39	63,70	69,97	79,60	88,65

X 2 « 1,992, X 3 = 0. Приближенное решение в предположении 3X 1 = X 2 , X 3 = 0 получено в работе [10]. Будем его называть случай 1, или случай кратных собственных чисел. Решение исходной задачи будем называть случай 2, или случаем некратных собственных чисел. Для случая 2 найти аналитическое решение не представляется возможным, поэтому получим его численно, используя решение случая 1 в качестве нулевого приближения.

Управление системой (4) – кусочно-постоянная функция с неизвестными временами переключения. Для всестороннего исследования пространства решений используется иерархическая стратегия оптимизации. Сначала проводится грубый поиск по сетке безразмерных времён переключения t1 , t2 продолжительности T, при этом параметры изменяются в пределах: t1 е[0,1; 1,0), 12 g(t1; 1,5), T g(12; 2,0], а шаг сетки составляет 0,02. В результате поиска и определяются возможные комбинации параметров, затем отбираются допустимые предварительные действительные решения системы (2), удовлетворяющие граничным условиям (4). В конце выполняется поиск решения в окрестности точек, полученных на предыдущем этапе. Новый поиск проводится по сетке с шагом 0,005 в пределах локального интервала ± 0,02 вокруг найденных ранее значений.

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

На рис. 1 представлены ф ( t ) , ф ( t ) для ряда допустимых амплитуд саккад. Для демонстрации влияния учета некратных корней на траектории системы рассмотрим изменения моментов переключений (долей времени постоянного управления (5)).

Из данных табл. 1 видно, что время этапа 1, затрачиваемое системой в случае 1, для движений с относительно малыми амплитудами было больше по сравнению со случаем 2. Для движений больших амплитуд это различие нивелируется с ростом амплитуды движения. На этапе 2 качественного различия между решениями не наблюдается. Этап 3, симметричный этапу 1, демон- стрировал противоположную тенденцию: система с кратными значениями затрачивала меньшее время на всех амплитудах, при этом наибольшее различие наблюдалось при наименьшей из рассмотренных амплитуд - 6,8°.

Таким образом, приближение для решения задачи (4), изложенное в статье [10], качественно соответствует полному решению задачи (4) и обладает близкими моментами переключения управления для движений с амплитудами более 15°. Для саккад меньших амплитуд расхождения решений значимы для ряда задач.

Традиционно для статистического описания используется главная последовательность: зависимость продолжительности саккады от ее амплитуды и максимальной скорости саккады от ее амплитуды. В литературе присутствует много вариантов математического описания главной последовательности.

Физиологическая основа ГП была первоначально изложена в [4; 5]. На этой базе были предложены различные модели с расширениями и вариантами параметров [22-28]. Однако из-за внутренней сложности явлений, которые описываются с помощью ГП, наблюдаются значительные различия в предлагаемых математических формулировках.

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

• I.    Продолжительность саккады от амплитуды

• 1.    Линейная модель:

T = T 0 A + T 1 , где коэффициенты имеют размерность T о , мс/°; T 1 , мс

• -    в статье [26]: до 12°: Т ₀ = 2 мс/°, T 1 = 20 мс;

• -    в работе [27]:    1,5 ° ~ 20 ° : Т ₀ = 2,9 мс/°,

• 2.    Квадратный корень:

Т 1 = 24 мс; 20 ° ~ 30 ° : Т ₀ = 1,3 мс/°, Т = 59 мс.

Т = Т о V A ; Т о , мс/7°

• -      в      источнике      [24]:       0 ° - 25 ° :

• 3.    Степенная модель: Т = Т₀ A m ; Т ₀, мс/ m °

Т _о е [ 11,8 мс/-^; 15,5 мс/7° ] .

• - в исследовании [27]:

• 1,5 ° - 20 ° : Т ₀ = 20 мс/ m °, m = 0,43; 20 ° - 30 ° :

Т ₀ = 40 мс/ m °, m = 0,24.

• II.    Соотношения для максимальной скорости и амплитуды

• 1.    Линейная модель: V P = V₀A + V 1 , где коэффициенты имеют размерность V о , c^-1; V 1 , °/c

• - в статье [26]: 5 ° ~ 25 ° : V ₀ = 185 c ^{- 1}, V 1 = 16,6 ° /c.

• 2.    Квадратичная зависимость V P = V ₀ JA ; V 0 , V° / c - в работе [27]:

1,5°- 20 ° : V ₀ = 100 V° /c;

20 °- 30 ° : V ₀ = 95 °F /c;

• - в источнике [28]:

• 3.    Степенная модель:

до 20°: V ₀ е[ 100 °°/c; 100,6 °°/c ] .

V_P = VAⁿ ; V ,

P 00 _n

• - в исследовании [27]:

• 1,5 ° - 20 ° :

• 4.    Экспоненциальная модель:

°/c                            °/c

V = 87 —, n = 0,56; 20° - 30°: V = 75 -^, 0n0n n = 0,57.

V _P = V 0 [ 1 - e ^- A / A ⁰ ] ; V 0 , ° /c, A ., °

• -    в статье [24]:

0 °- 25 ° :

V ₀ е [ 557,5 °/c; 614,19 °/c ] , A ₀ е [ 5,6 ° ; 9,3 ° ] ;

• -    в работе [27]: 1,5 ° - 20 ° : V ₀ = 442 ° /c,

A ₀ = 7,1 ° ; 20°- 30 ° : V ₀ = 759 ° /c, A ₀ = 25,6 ° .

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

Поскольку получение явных аналитических соотношений, основанное на решении задачи (4), не представляется возможным, в последующем анализе используются численные решения. Для анализа сформировано 34 модельных траектории с амплитудами A i в диапазоне от 0° до 40°: A i е [ 0°; 40° ] , i = 1...34. Эффективность каждой из моделей ГП оценивалась по среднему абсолютному отклонению ( mean absolute error , MAE ), также будем его назвать невязкой. MAE отражает среднюю величину абсолютной ошибки между предсказанными моделью и фактическими значениями параметров. Критерий точности имеет следующую формулировку:

MAE = 1 X П = 1 Т ( A ) - ^Т ( A i ) ,

Таблица 3

Коэффициенты функциональных соотношений

ГП, вычисленные по модельным данным

Линейная	V 0	193,24	T 0	1,06
	V 1	18,47	T 1	48,34
	MAE ( V P )	4,62	MAE ( T )	1,20
Степенная	V 0	77,92	T 0	26,61
	n	0,66	m	0,33
	MAE ( V P )	0,53	MAE ( T )	0,96
Экспоненциальная	V 0	716,99	–	–
	A 0	12,46
	MAE ( V P )	34,10

1000 ■

800 ■

600 ■

400 ■ о о. о к я и о

Данные модели-некратных

Данные модели-кратных Степенная функция-некратных

Линейная-некратных

Экспоненциальная-некратных

20 25 30 35 40

Амплитуда саккад,0

Рис. 2. Главная последовательность пиковых скоро-стей-амплитуды и продолжительности-амплитуды

MAE- = 1 Z ”„ V- (Ai)-VP (-A')| • где n - количество точек; T (Ai), VP (A) - показатели, рассчитанные для Ai в соответствии с решением задачи (4), T (A), Vp (A) - i-е предсказанные значения для различных моделей основной последовательности.

Для каждой модели предыдущего пункта был про- веден поиск коэффициентов моделей ГП. Результаты приведены в табл. 3 и на рис. 2.

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

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

Для получения максимального отклонения в вычислении максимальной угловой скорости выбиралось равномерное разбиение по времени со сдвигом относительно начала модельной траектории:

Ow = Teal — dt /2) — П • dt, где Tideal – момент времени, соответствующий пиковой скорости в непрерывном времени, dt = 1/f – шаг дискретизации, n – целое число, чтобы вычислить первую точку разбиения (отступ от нуля). Такое разбиение по времени гарантированно дает максимальную ошибку при расчете максимальной угловой скорости глаза на саккаде. Обозначим значения MAE для каждой модели главной последовательности при разных частотах дискретизации как MAE(T, f), MAE(VP, f), f – частота дискретизации.

Согласно табл. 4, линейная модель имеет наименьшую ошибку при моделировании данных с частотой дискретизации 30 Гц. Стоит учитывать, что на саккаду приходится до четырех точек, причем для саккад до ~15° приходится лишь три точки. В этом случае применение, например, стандартных видов сглаживания значимо искажает полученные оценки, и классификация точки как принадлежащей к саккаде или нет имеет низкую степень достоверности. Таким образом, ГП, полученная на данных с низкой величиной дискретизации, не может считаться близкой к истинной. С увеличением частоты наилучшее описание – степенная зависимость, как для соотношения «максимальная скорость – амплитуда», так и для соотношения «продолжительность – амплитуда». Экспоненциальная зависимость не приводится в табл. 4, так как демонстрирует невязки порядка 40 мс для продолжительности саккад и 40°/c для скоростей для всех рассмотренных частот дискретизации.

Таблица 4

Сравнение математических описаний ГП на модельных данных при различной частоте дискретизации

Модели оценки максимальной скорости
	Линейная			Степенная
Hz		V 1	MAE(Vp,f)		п	MAE(V ,j
30	0,10	15,54	13,10	4,05	1,43	38,53
60	93,25	18,13	16,12	45,68	0,78	6,72
90	122,62	18,18	16,16	60,40	0,71	2,54
100	127,07	18,19	16,17	62,74	0,70	2,01
250	143,01	18,22	16,19	71,27	0,67	0,87
500	138,24	18,50	18,20	72,59	0,66	1,02
1000	152,25	17,98	14,64	73,04	0,66	1,15
1250	145,92	18,22	16,20	72,86	0,66	1,03
Модели оценки продолжительности саккады
Линейная				Степенная
Hz		т 1	MAE(jj)	Го	m	MAE(TJ)
30	37,56	0,93	1,03	27,10	0,25	0,41
60	37,38	1,45	4,73	22,12	0,38	4,65
90	43,07	1,27	3,52	27,29	0,32	2,99
100	42,36	1,33	4,16	25,30	0,35	3,28
250	46,12	1,32	3,24	28,69	0,32	2,16
500	47,36	1,30	3,06	30,39	0,31	2,08
1000	48,03	1,28	3,09	30,73	0,30	2,15
1250	48,05	1,29	3,08	30,71	0,31	2,09

Рис. 3. Главная последовательность «максимальная ско- рость – амплитуда» и «длительность – амплитуда»

Заключение

Итак, в работе получено математическое описание саккадического движения глаза в различных аспектах.

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

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

T = T ₀ A^m ,

V P = V 0 Aⁿ .

При изменении частоты от 90 до 1250 Гц для соотношения «максимальная скорость – амплитуда» коэффициент масштаба ( V 0 ) принадлежит множеству [0,85; 2,54], а степень ( n ) принадлежит множеству [0,66; 0,71]; для соотношения «продолжительность – амплитуда»

коэффициент масштаба ( T 0 ) принадлежит множеству [27; 31], а степень ( m ) принадлежит множеству [0,31; 0,32]. Полученные значения согласуются с приводимыми в ранее опубликованных исследованиях.

Прикладной результат исследования заключается в том, что при описании экспериментальных данных саккадической главной последовательностью необходимо указывать частоту съема данных, при необходимости прямое сравнение корректно проводить с результатами, полученными или приведенными к той же частоте. Для сравнения различных исследований представляется применимым подход сравнения с параметрами модельных саккад (например, приведенными в табл. 4 и рис. 3 данной работы) или параметрами, полученными аналогичным образом на основе математической модели, например, (4) с приводимыми в данной работе параметрами. Стоит отметить, что точность моделей зависит от принимаемых допущений и для амплитуд саккад более 10–15° результаты упрощённой модели практически совпадают с результатами полной модели. При исследовании микросаккад стоит провести дополнительную оценку точности используемой модели и отдать предпочтение моделями типа (4).