Журнал вычислительной математики и математической физики, 2021, T. 61, № 7, стр. 1162-1171

Распознавание квазипериодической последовательности, включающей неизвестное число нелинейно-растянутых эталонных подпоследовательностей

А. В. Кельманов^1, 2, Л. В. Михайлова^1, *, П. С. Рузанкин^1, 2, **, С. А. Хамидуллин¹

¹ Ин-т матем. им. С.Л. Соболева
630090 Новосибирск, пр-т акад. Коптюга, 4, Россия

² Новосибирский гос. ун-т
630090 Новосибирск, ул. Пирогова, 2, Россия

^* E-mail: mikh@math.nsc.ru
^** E-mail: ruzankin@math.nsc.ru

Поступила в редакцию 26.11.2020
После доработки 26.11.2020
Принята к публикации 11.03.2021

DOI: 10.31857/S0044466921070097

Полный текст (PDF)

Аннотация

Рассматривается неизученная экстремальная задача, которая индуцируется одной из задач помехоустойчивого распознавания квазипериодической последовательности, а именно, задачей распознавания последовательности $Y$ длины $N$ как последовательности, порожденной некоторой последовательностью $U$, принадлежащей заданному конечному множеству $W$ (алфавиту) последовательностей. Каждая последовательность $U$ из $W$ порождает экспоненциальное по мощности множество $\mathcal{X}(U)$ последовательностей, объединяющее все последовательности длины $N$, которые в качестве подпоследовательностей включают переменное число допустимых квазипериодических (флуктуационных) повторов последовательности $U$. Каждый квазипериодический повтор порождается допустимыми преобразованиями последовательности $U$, а именно, сдвигами и растяжениями. Задача распознавания состоит в выборе последовательности $U$ из $W$ и аппроксимации последовательности $Y$ элементом $X$ из множества $\mathcal{X}(U)$ последовательностей. Критерием аппроксимации является минимум суммы квадратов расстояний между элементами последовательностей. Мы показываем, что рассматриваемая задача эквивалентна задаче суммирования элементов двух числовых последовательностей, в которой требуется минимизировать сумму неизвестного числа $M$ слагаемых, каждое из которых является разностью невзвешенной автосвертки растянутой на переменную длину последовательности $U$ (путем кратных повторов ее элементов) и взвешенной свертки этой растянутой последовательности с подпоследовательностью из $Y$. Мы доказываем, что рассматриваемая экстремальная задача и вместе с ней задача распознавания разрешимы за полиномиальное время. Примерами численного моделирования проиллюстрирована применимость алгоритма к решению модельных прикладных задач помехоустойчивой обработки ECG-подобных и PPG-подобных квазипериодических сигналов (electrocardiogram-like and photoplethysmogram-like signals). Библ. 9. Фиг. 5.

Ключевые слова: числовые последовательности, распознавание, квазипериодическая последовательность, полиномиальная разрешимость, разность взвешенных сверток.

ВВЕДЕНИЕ

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

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

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

1. ФОРМУЛИРОВКА ЗАДАЧИ, ЕЕ ИСТОКИ И ТРАКТОВКИ

Рассматриваемая экстремальная задача имеет следующую формулировку.

Задача 1. Дано: числовая последовательность $Y = ({{y}_{1}}, \ldots ,{{y}_{N}})$, совокупность $W = \{ {{U}^{{(1)}}}, \ldots ,{{U}^{{(K)}}}\,{\text{|}}\,{{U}^{{(k)}}} = (u_{1}^{{(k)}}, \ldots ,u_{{{{q}_{k}}}}^{{(k)}}),k = 1, \ldots ,K\} $ $W = \{ {{U}^{{(1)}}}, \ldots ,{{U}^{{(K)}}}\,{\text{|}}\,{{U}^{{(k)}}} = (u_{1}^{{(k)}}, \ldots ,u_{{{{q}_{k}}}}^{{(k)}}),k = 1, \ldots ,K\} $, натуральные числа ${{T}_{{{\text{max}}}}}$ и $\ell $. Найти: числовую последовательность $U = ({{u}_{1}}, \ldots ,{{u}_{{q(U)}}}) \in W$, набор $\mathcal{M} = \{ {{n}_{1}}, \ldots {{n}_{m}}, \ldots \} $ номеров последовательности $Y$, набор $\mathcal{P} = \{ {{p}_{1}}, \ldots ,{{p}_{m}}, \ldots \} $ натуральных чисел, набор $\mathcal{J} = \{ {{J}^{{(1)}}}, \ldots ,{{J}^{{(m)}}}, \ldots \} $ сжимающих отображений, в котором ${{J}^{{(m)}}}:\{ 1, \ldots ,{{p}_{m}}\} \to \{ 1, \ldots ,q(U)\} $, а также размерность $M$ этих наборов, доставляющих минимум целевой функции

(1.1)

$F(\mathcal{U},\mathcal{M},\mathcal{P},\mathcal{J}) = \sum\limits_{m = 1}^M \,\sum\limits_{i = 1}^{{{p}_{m}}} \,\{ u_{{{{J}^{{(m)}}}(i)}}^{2} - 2{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}\} ,$

при ограничениях

(1.2)

$q(U) \leqslant {{p}_{m}} \leqslant \ell \leqslant {{T}_{{{\text{max}}}}} \leqslant N,\quad m = 1, \ldots ,M,$

(1.3)

${{p}_{{m - 1}}} \leqslant {{n}_{m}} - {{n}_{{m - 1}}} \leqslant {{T}_{{{\text{max}}}}},\quad m = 2, \ldots ,M,$

(1.4)

${{p}_{M}} \leqslant N - {{n}_{M}} + 1,$

на элементы искомых наборов $\mathcal{M}$, $\mathcal{P}$, и при ограничениях

(1.5)

$\begin{array}{*{20}{c}} {{{J}^{{(m)}}}(1) = 1,\quad {{J}^{{(m)}}}({{p}_{m}}) = q(U),} \\ {0 \leqslant {{J}^{{(m)}}}(i) - {{J}^{{(m)}}}(i - 1) \leqslant 1,\quad i = 2, \ldots ,{{p}_{m}},\quad m = 1, \ldots ,M,} \end{array}$

на элементы искомых сжимающих отображений.

Приведем несколько трактовок задачи 1.

Из формулировки задачи 1 и вида целевой функции (1.1) следует, что задача 1 – задача об оптимальном (в смысле минимума (1.1)) суммировании элементов двух последовательностей, одна из которых $Y$ задана, а другая $U$ принадлежит заданному множеству последовательностей. Эквивалентная перезапись целевой функции (1.1) в виде

$F(\mathcal{U},\mathcal{M},\mathcal{P},\mathcal{J}) = \sum\limits_{m = 1}^M \,\left\{ {\sum\limits_{i = 1}^{{{p}_{m}}} \,u_{{{{J}^{{(m)}}}(i)}}^{2} - 2\sum\limits_{i = 1}^{{{p}_{m}}} \,{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}} \right\}$

позволяет трактовать задачу 1 как задачу минимизации суммы разностей взвешенных сверток. Действительно, при каждом $m = 1, \ldots ,M$ выражение $\sum\nolimits_{i = 1}^{{{p}_{m}}} \,u_{{{{J}^{{(m)}}}(i)}}^{2}$ в фигурных скобках – невзвешенная автосвертка последовательности ${{u}_{{{{J}^{{(m)}}}(i)}}}$, $i = 1, \ldots ,{{p}_{m}}$, полученной из элементов $U$ растяжением путем дублирования элементов, a выражение $\sum\nolimits_{i = 1}^{{{p}_{m}}} \,{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}$ – свертка этой растянутой последовательности с подпоследовательностью из $Y$, имеющей ту же длину ${{p}_{m}}$ (коэффициент 2 – вес этой свертки).

Другая возможная трактовка – задача совместного выбора последовательности $U \in W$ и аппроксимации последовательности $Y$ последовательностью $X \in \mathcal{X}(U)$, по критерию минимума суммы квадратов расстояний между элементами $Y$ и $X$, т.е. задача

(1.6)

${{\left\| {Y - X} \right\|}^{2}} \to \mathop {min}\limits_{U,\mathcal{X}(U)} .$

Здесь $\mathcal{X}(U)$, $U \in W$, – множество допустимых аппроксимирующих последовательностей, порожденных последовательностью $U$. Каждый элемент $X = ({{x}_{1}}, \ldots ,{{x}_{N}}) \in \mathcal{X}(U)$ однозначно определяется наборами $\mathcal{M}$, $\mathcal{P}$, $\mathcal{J}$, удовлетворяющими ограничениям (1.2)–(1.5), по правилу

(1.7)

${{x}_{n}} = \sum\limits_{m = 1}^M \,h_{{n - {{n}_{m}} + 1}}^{{(m)}},\quad n = 1, \ldots ,N,$

где

(1.8)

$h_{i}^{{(m)}} = \left\{ \begin{gathered} {{u}_{{{{J}^{{(m)}}}(i)}}},\quad {\text{если}}\quad i = 1, \ldots ,{{p}_{m}}, \hfill \\ 0,\quad {\text{если}}\quad i < 1,\quad i > {{p}_{m}},\quad m = 1, \ldots ,M. \hfill \\ \end{gathered} \right.$

Равенство (1.8) устанавливает связь между элементами последовательностей $h_{i}^{{(m)}}$ и $U$. Из этого равенства видно, что каждая из последовательностей $h_{i}^{{(m)}}$ – растянутая до длины ${{p}_{m}}$ последовательность $U$, а кратность дублирования ее элементов определяется формулой

(1.9)

$k_{t}^{{(m)}} = {\text{|}}\{ i\,|\,{{J}^{{(m)}}}(i) = t,i \in \{ 1, \ldots ,{{p}_{m}}\} \} {\text{|}},\quad t = 1, \ldots ,q(U),$

причем

(1.10)

${{p}_{m}} = k_{1}^{{(m)}} + \ldots + k_{{q(U)}}^{{(m)}},\quad m = 1, \ldots ,M.$

Формула (1.7) – сумма $M$ растянутых последовательностей вида (1.8). Таким образом, последовательность $X$ включает в себя $M$ повторов растянутых последовательностей $U$. Значение индекса $n = {{n}_{m}}$, ${{n}_{m}} \in \mathcal{M}$, определяет начальный номер $m$-го повтора, а значение $p = {{p}_{m}}$, ${{p}_{m}} \in \mathcal{P}$, и отображение $J = {{J}^{{(m)}}}$, ${{J}^{{(m)}}} \in \mathcal{J}$, – длину этого повтора и кратности дублирования элементов из $U$ в этом повторе.

Неформально можно сказать, что каждая последовательность $X \in \mathcal{X}(U)$ включает в себя неизвестное число допустимых квазипериодических повторов последовательности $U$. При этом каждый повтор определяется 1) сдвигом $U$ на переменную величину, которая между соседними повторами не превышает ${{T}_{{{\text{max}}}}} \leqslant N$; 2) допустимым растяжением последовательности $U$, путем дублирования ее компонент.

Из (1.7) и (1.8) следует, что для каждого $U \in W$ между элементами множества $\mathcal{X}(U)$ и наборами $\mathcal{M}$, $\mathcal{P}$, $\mathcal{J}$, удовлетворяющими ограничениям (1.2)–(1.5), существует взаимооднозначное соответствие, т.е. $X = X(U,\mathcal{M},\mathcal{P},\mathcal{J})$, поэтому задачу (1.6) можно эквивалентно записать в виде

(1.11)

${{\left\| {Y - X} \right\|}^{2}} = {{\left\| {Y - X(U,\mathcal{M},\mathcal{P},\mathcal{J})} \right\|}^{2}} \to \mathop {min}\limits_{U,\mathcal{M},\mathcal{P},\mathcal{J}} .$

Наконец, преобразовав ${{\left\| {Y - X} \right\|}^{2}}$ с учетом (1.7) и (1.8), имеем

$\sum\limits_{n = 1}^N \,{{({{x}_{n}} - {{y}_{n}})}^{2}} = \sum\limits_{n = 1}^N \,y_{n}^{2} + \sum\limits_{m = 1}^M \,\sum\limits_{i = 1}^{{{p}_{m}}} \,\{ u_{{{{J}^{{(m)}}}(i)}}^{2} - 2{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}\} .$

Легко видеть, что первое слагаемое в правой части этого равенства – константа, не зависящая от переменных задачи 1, а второе совпадает с целевой функцией задачи 1. Поэтому задача (1.11), а вместе с ней и задача (1.6), эквивалентна оптимизационной задаче 1.

Чтобы проиллюстрировать данную интерпретацию, на фиг. 1 и 2 в виде графиков изображены входные данные задачи 1. Здесь и далее считаем, что номерам элементов последовательностей из алфавита и последовательности $Y$ соответствуют дискретно-временные отсчеты непрерывных сигналов. На фиг. 1 приведен пример алфавита $W$, содержащего 6 последовательностей, на фиг. 2 – пример последовательности $Y$, порожденной одной из последовательностей алфавита.

Фиг. 1.

Пример входных значений задачи 1. Алфавит $W$.

Фиг. 2.

Пример входных значений задачи 1. Последовательность $Y$, подлежащая обработке.

На фиг. 3a изображена последовательность $U$ из алфавита $W$, на фиг. 3б – пример последовательности $X \in \mathcal{X}(U)$, построенной по правилам (1.7), (1.8), (1.9), (1.10) для некоторых допустимых значений $\{ {{p}_{m}}\} $ и $\{ k_{t}^{{(m)}}\} $. Ступенчатость сигнала на этой фигуре обусловлена кратными повторениями элементов из $U$.

Фиг. 3.

Пример последовательности $U$ из алфавита $W$ и одной из последовательностей $X \in \mathcal{X}(U)$.

2. БЛИЗКИЕ ПО ПОСТАНОВКЕ ЗАДАЧИ

Задача 1 является обобщением ранее исследованных задач распознавания. Частный случай задачи 1, в котором ${{q}_{1}} = \ldots = {{q}_{K}} = q$, ${{p}_{m}} = q$, ${{J}^{{(m)}}}(i) = i$, $m = 1, \ldots ,M$, исследовался в [1]. В [2] была рассмотрена модификация этого частного случая, когда $M$ является частью входа задачи. В этом частном случае и его модификации, как и в задаче 1, требуется распознать квазипериодическую последовательность, но все повторы порождающей последовательности в ней идентичны. В [1] и [2] были построены алгоритмы, позволяющие получить оптимальное решение за время $\mathcal{O}(K{{T}_{{{\text{max}}}}}N)$ и $\mathcal{O}(KM{{T}_{{{\text{max}}}}}N)$. Другой частный случай, когда $K = 1$, был исследован в [3]. В этом частном случае алфавит состоит из единственной последовательности, поэтому задачу можно интерпретировать как задачу аппроксимации. В этой работе приведен точный алгоритм, позволяющий получить точное решение за время $\mathcal{O}(T_{{{\text{max}}}}^{3}N)$.

Близкими в постановочном плане являются задачи распознавания квазипериодических последовательностей по усеченным данным [4] и [5]. В этих задачах, как и в задаче 1, предполагалось, что распознаваемая квазипериодическая последовательность включает в себя искаженные допустимым образом повторы некоторой последовательности, принадлежащей заданному множеству. В отличие от задачи 1, здесь предполагалось, что множество допустимых преобразований последовательности состоит из всевозможных ее подпоследовательностей, полученных удалением некоторого количества начальных и конечных отсчетов. Для этих задач в цитируемых работах получены точные полиномиальные алгоритмы для их решения.

Алгоритм решения задачи 1 является подходящим инструментом для решения прикладных задач распознавания и анализа сигналов, имеющих квазипериодическую структуру в виде флуктуирующих повторов участков сигнала. Распознавание и анализ таких сигналов актуальны для различных приложений, имеющих дело с обработкой импульсных сигналов, полученных от природных источников: биомедицинские, геофизические и т.п. Далее, в разделе численное моделирование, будут приведены примеры обработки медицинских сигналов, имеющих вид квазипериодических последовательностей (ECG, PPG).

3. РАЗМЕР МНОЖЕСТВА ДОПУСТИМЫХ РЕШЕНИЙ

Приведенная выше интерпретация позволяет заметить, что число допустимых решений задачи 1 совпадает с мощностью множества $\mathcal{X} = \bigcup\nolimits_{U \in W}^{} {\mathcal{X}(U)} $. Из комбинаторных соображений легко видеть, что при каждом $U \in W$ справедлива оценка

$\begin{gathered} \left| {\mathcal{X}(U)} \right| \leqslant (N - q(U) + 1)\sum\limits_{M = 1}^{{{M}_{{{\text{max}}}}}(U)} \,q{{(U)}^{{({{T}_{{{\text{max}}}}} - q(U))M}}}{{({{T}_{{{\text{max}}}}} - q(U) + 1)}^{{2M - 1}}} \leqslant \\ \leqslant q{{(U)}^{{({{T}_{{{\text{max}}}}} - q(U)){{M}_{{{\text{max}}}}}(U)}}}(N - q(U) + 1){{({{T}_{{{\text{max}}}}} - q(U) + 1)}^{{2{{M}_{{{\text{max}}}}}(U) - 1}}}{{M}_{{{\text{max}}}}}(U), \\ \end{gathered} $

где ${{M}_{{{\text{max}}}}}(U) = \left\lfloor {N{\text{/}}q(U)} \right\rfloor $ – максимально возможное число повторов последовательности $U$ в $Y$. Данная верхняя оценка позволяет оценить влияние, которое параметры задачи 1 оказывают на мощность множества допустимых решений.

С другой стороны, за исключением тривиального случая ${{T}_{{{\text{max}}}}} = q(U)$ справедлива нижняя оценка

$\left| {\mathcal{X}(U)} \right| \geqslant {{2}^{{\left\lfloor {\tfrac{{N - q(U) + 1}}{{q(U) + 1}}} \right\rfloor }}},\quad U \in W.$

Отсюда следует, что

$\left| \mathcal{X} \right| = \sum\limits_{U \in W} \,\left| {\mathcal{X}(U)} \right| \geqslant K{{2}^{{\left\lfloor {\tfrac{{N - {{q}_{{{\text{max}}}}} + 1}}{{{{q}_{{{\text{max}}}}} + 1}}} \right\rfloor }}},$

где ${{q}_{{{\text{max}}}}} = \mathop {max}\limits_{U \in W} q(U)$, а $K$ – мощность алфавита последовательностей. Это означает, что если ${{q}_{{{\text{max}}}}}$ ограничено некоторой константой (что типично для приложений), мощность множества $\mathcal{X}$ допустимых решений задачи 1 экспоненциально растет с ростом $N$. Очевидно, что перебрать напрямую элементы этого множества за приемлемое время вряд ли удастся, так как $N$ – часть входа задачи. Несмотря на этот экспоненциальный рост, ниже будет приведен алгоритм, позволяющий получить оптимальное решение за полиномиальное время.

4. ОСНОВЫ АЛГОРИТМА

Для построения алгоритма решения задачи 1 нам потребуется следующая вспомогательная

Задача 2. Дано: числовые последовательности $Y = ({{y}_{1}}, \ldots ,{{y}_{N}})$, $U = ({{u}_{1}}, \ldots ,{{u}_{{q(U)}}})$, и натуральные числа ${{T}_{{{\text{max}}}}}$, $\ell $. Найти: набор $\mathcal{M} = \{ {{n}_{1}}, \ldots {{n}_{m}}, \ldots \} $ номеров последовательности $Y$, набор $\mathcal{P} = \{ {{p}_{1}}, \ldots ,{{p}_{m}}, \ldots \} $ натуральных чисел, набор $\mathcal{J} = \{ {{J}^{{(1)}}}, \ldots ,{{J}^{{(m)}}}, \ldots \} $ сжимающих отображений, в котором ${{J}^{{(m)}}}:\{ 1, \ldots ,{{p}_{m}}\} \to \{ 1, \ldots ,q(U)\} $, а также размерность $M$ этих наборов, которые минимизируют целевую функцию

(4.1)

$G(\mathcal{M},\mathcal{P},\mathcal{J}) = F( \bullet \,{\text{|}}\,U)$

при ограничениях (1.2), (1.3), (1.4) на элементы искомых наборов $\mathcal{M}$, $\mathcal{P}$, и при ограничениях (1.5) на элементы искомых сжимающих отображений.

Для полноты изложения приведем алгоритм ее решения, базирующийся на результатах, полученных в [3]. Проанализировав ограничения (1.2), (1.3), (1.4), входящие в условия задачи 2, легко видеть, что справедливо

Утверждение 1. Пусть выполнены условия задачи 2, тогда для компонент наборов $\mathcal{M}$ и $\mathcal{P}$ справедливо:

1) ${{n}_{m}} \in \omega $, $m = 1, \ldots ,M$, где

$\omega = \{ 1, \ldots ,N - q(U) + 1\} ;$

2) ${{n}_{m}} \in {{\omega }^{ + }}$, $m = 2, \ldots ,M$, где

${{\omega }^{ + }} = \{ q(U) + 1, \ldots ,N - q(U) + 1\} ;$

3) если ${{n}_{m}} = n$, $m = 1, \ldots M$, $n \in \omega $, то ${{p}_{m}} \in \delta (n)$, где

$\delta (n) = \{ q(U), \ldots ,\min \{ \ell ,N - n + 1\} \} ;$

4) если ${{n}_{m}} = n$, $m = 2, \ldots ,M$, $n \in {{\omega }^{ + }}$, то ${{n}_{{m - 1}}} \in \gamma (n)$, где

$\gamma (n) = \{ k\,{\text{|}}\,\max \{ n - {{T}_{{\max }}},1\} \leqslant k \leqslant n - q(U)\} ;$

5) если ${{n}_{m}} = n$ и ${{n}_{{m - 1}}} = j$, $m = 2, \ldots ,M$, $n \in {{\omega }^{ + }}$, $j \in \gamma (n)$ то ${{p}_{{m - 1}}} \in \theta (n,j)$, где

$\theta (n,j) = \{ q(U), \ldots ,\min \{ \ell ,n - j\} \} .$

Опираясь на это утверждение, выпишем алгоритм решения задачи 2.

Алгоритм ${{\mathcal{A}}_{1}}$

Вход: $Y$, $U$, ${{T}_{{{\text{max}}}}}$ и $\ell $.

Прямой ход алгоритма.

Шаг 1. Для каждого $n = 1, \ldots ,N - q(U) + 1$ выполним:

Для каждого $p = q(U), \ldots ,\min \{ \ell ,N - n + 1\} $ выполним:

(1) Вычислим

${{w}_{{s,t}}} = u_{t}^{2} - 2{{y}_{{n + s - 1}}}{{u}_{t}},\quad s = 1, \ldots ,p,\quad t = 1, \ldots ,q(U).$

(2) Найдем значение $W{\kern 1pt} *$ по формулам

$W{\kern 1pt} * = {{W}_{{p,q(U)}}},$

${{W}_{{s,t}}} = \min \{ {{W}_{{s - 1,t}}},{{W}_{{s - 1,t - 1}}}\} + {{w}_{{s,t}}},\quad s = 1, \ldots ,p,\quad t = 1, \ldots ,q(U),$

${{W}_{{s,t}}} = \left\{ \begin{gathered} 0,\quad s = 0,\quad t = 0, \hfill \\ + \infty ,\quad s = 0,\quad t = 1, \ldots ,q(U), \hfill \\ + \infty ,\quad s = 1, \ldots ,p,\quad t = 0, \hfill \\ \end{gathered} \right.$

и последовательность $J{\kern 1pt} *{\kern 1pt} (1), \ldots ,J{\kern 1pt} *{\kern 1pt} (p)$ по формулам

$\begin{array}{*{20}{c}} {J{\kern 1pt} *{\kern 1pt} (p) = q(U);} \\ {J{\kern 1pt} *{\kern 1pt} (i - 1) = \left\{ \begin{gathered} J{\kern 1pt} *{\kern 1pt} (i),\quad {\text{если}}\quad {{W}_{{i - 1,J{\kern 1pt} *(i)}}} \leqslant {{W}_{{i - 1,J{\kern 1pt} *(i) - 1}}}, \hfill \\ J{\kern 1pt} *{\kern 1pt} (i) - 1,\quad {\text{если}}\quad {{W}_{{i - 1,J{\kern 1pt} *(i)}}} > {{W}_{{i - 1,J{\kern 1pt} *(i) - 1}}}, \hfill \\ \end{gathered} \right.} \\ {i = p,p - 1, \ldots ,2.} \end{array}$

(3) Положим $W{\kern 1pt} *{\kern 1pt} (n,p) = W{\kern 1pt} *$; $J{\kern 1pt} *{\kern 1pt} (n,p) = \{ J{\kern 1pt} *{\kern 1pt} (i),i = 1, \ldots ,p\} $.

Шаг 2. Вычислим совокупность значений $G(n,p)$, $p = q(U), \ldots ,\ell $, $n = 1, \ldots ,N - p + 1$, по формулам

$G(n,p) = \left\{ \begin{gathered} W{\kern 1pt} *{\kern 1pt} (n,p),\quad n \in \omega ,\quad p \in \delta (n), \hfill \\ \min \{ 0,\mathop {\min }\limits_{j \in \gamma (n)} \mathop {\min }\limits_{i \in \theta (n,j)} G(j,i)\} + W{\kern 1pt} *{\kern 1pt} (n,p)),\quad n \in {{\omega }^{ + }},\quad p \in \delta (n), \hfill \\ \end{gathered} \right.$

а так же значение $G{\kern 1pt} *$ по формуле

$G{\kern 1pt} * = \mathop {min}\limits_{n \in \omega } \mathop {min}\limits_{p \in \delta (n)} G(n,p).$

Положим ${{F}_{A}} = G{\kern 1pt} *$.

Обратный ход алгоритма.

Шаг 3. Вычислим $\pi (n)$ и $I(n)$, $n = 1, \ldots ,N - q(U) + 1$, по формулам

$I(n) = \left\{ \begin{gathered} 0,\quad {\text{если}}\quad n \in \omega {{\backslash }}{{\omega }^{ + }}, \hfill \\ 0,\quad {\text{если}}\quad \mathop {min}\limits_{j \in \gamma (n)} \mathop {min}\limits_{i \in \theta (n,j)} G(j,i) \geqslant 0,\quad n \in {{\omega }^{ + }}, \hfill \\ arg\mathop {min}\limits_{j \in \gamma (n)} \left\{ {\mathop {min}\limits_{i \in \theta (n,j)} G(j,i)} \right\},\quad {\text{если}}\quad \mathop {min}\limits_{j \in \gamma (n)} \mathop {min}\limits_{i \in \theta (n,j)} G(j,i) < 0,\quad n \in {{\omega }^{ + }}, \hfill \\ \end{gathered} \right.$

$\pi (n) = \left\{ \begin{gathered} 0,\quad {\text{если}}\quad I(n) = 0, \hfill \\ arg\mathop {min}\limits_{i \in \theta (n,I(n))} G(I(n),i),\quad {\text{если}}\quad I(n) > 0,\quad n \in \omega {\kern 1pt} . \hfill \\ \end{gathered} \right.$

Шаг 4. Найдем компоненты вспомогательных наборов $({{\pi }_{1}}, \ldots ,{{\pi }_{M}})$ и $({{\nu }_{1}}, \ldots ,{{\nu }_{M}})$ и их размерность $M$ по формулам

${{\pi }_{1}} = \arg \mathop {\min }\limits_{p \in \{ q(U), \ldots ,\ell \} } \left\{ {\mathop {\min }\limits_{n \in \{ 1, \ldots ,N - p + 1\} } G(n,p)} \right\},$

${{\nu }_{1}} = arg\mathop {min}\limits_{n \in \{ 1, \ldots ,N - {{\pi }_{1}} + 1\} } G(n,{{\pi }_{1}}),$

${{\pi }_{m}} = \pi ({{\nu }_{{m - 1}}}),\quad {{\nu }_{m}} = I({{\nu }_{{m - 1}}}),\quad m = 2, \ldots ,M,$

где $M$ – наименьшее значение $m$, при котором $\pi ({{\nu }_{m}}) = 0$.

Шаг 5. Положим: ${{M}_{{{{A}_{1}}}}} = M$, ${{\mathcal{M}}_{{{{A}_{1}}}}} = \{ {{\nu }_{M}}, \ldots ,{{\nu }_{1}}\} $, ${{\mathcal{P}}_{{{{A}_{1}}}}} = \{ {{\pi }_{M}}, \ldots ,{{\pi }_{1}}\} $, ${{\mathcal{J}}_{{{{A}_{1}}}}} = \{ J{\kern 1pt} *{\kern 1pt} ({{\nu }_{M}},{{\pi }_{M}}), \ldots ,J{\kern 1pt} *{\kern 1pt} ({{\nu }_{1}},{{\pi }_{1}})\} .$

Выход: ${{M}_{{{{A}_{1}}}}}$, ${{\mathcal{M}}_{{{{A}_{1}}}}}$, ${{\mathcal{P}}_{{{{A}_{1}}}}}$, ${{\mathcal{J}}_{{{{A}_{1}}}}}$ и ${{G}_{{{{A}_{1}}}}}$.

Утверждение 2 (см. [3]). Алгоритм ${{\mathcal{A}}_{1}}$ находит точное решение задачи 2 за время $\mathcal{O}(T_{{{\text{max}}}}^{3}N)$.

5. АЛГОРИТМ

Опираясь на вспомогательную задачу 2 и утверждение 2, сформулируем алгоритм решения задачи 1.

Алгоритм $\mathcal{A}$

Вход: $Y$, $W$, ${{T}_{{{\text{max}}}}}$ и $\ell $.

Шаг 1. Для каждого $U \in W$ выполним алгоритм ${{\mathcal{A}}_{1}}$ с входными данными $Y$, $U$, ${{T}_{{{\text{max}}}}}$ и $\ell $. Положим $M(U) = {{M}_{{{{A}_{1}}}}}$, $\mathcal{M}(U) = {{\mathcal{M}}_{{{{A}_{1}}}}}$, $\mathcal{P}(U) = {{\mathcal{P}}_{{{{A}_{1}}}}}$, $\mathcal{J}(U) = {{\mathcal{J}}_{{{{A}_{1}}}}}$, ${{G}_{U}} = {{G}_{{{{A}_{1}}}}}$.

Шаг 2. Вычислим

(5.1)

$\begin{gathered} {{F}_{A}} = \mathop {min}\limits_{U \in W} {{G}_{U}}, \\ {{U}_{A}} = arg\mathop {min}\limits_{U \in W} {{G}_{U}}. \\ \end{gathered} $

Шаг 3. Положим ${{M}_{A}} = M({{U}_{A}})$, ${{\mathcal{M}}_{A}} = \mathcal{M}({{U}_{A}})$, ${{\mathcal{P}}_{A}} = \mathcal{P}({{U}_{A}})$, ${{\mathcal{J}}_{A}} = \mathcal{J}({{U}_{A}})$.

Выход: ${{U}_{A}}$, ${{M}_{A}}$, ${{\mathcal{M}}_{A}}$, ${{\mathcal{P}}_{A}}$, ${{\mathcal{J}}_{A}}$ и ${{F}_{A}}$.

Теорема 1. Алгоритм $\mathcal{A}$ находит точное решение задачи 1 за время $\mathcal{O}(KT_{{{\text{max}}}}^{3}N)$.

Доказательство. Оптимальность решения, найденного алгоритмом $\mathcal{A}$, следует из Утверждения 2 и следующей цепочки равенств:

$\begin{gathered} F{\kern 1pt} *\;\mathop = \limits^{\text{1}} \;\mathop {min}\limits_{U,\mathcal{M},\mathcal{P},\mathcal{J}} F(U,\mathcal{M},\mathcal{P},\mathcal{J})\;\mathop = \limits^{\text{2}} \;\mathop {min}\limits_U \mathop {min}\limits_{\mathcal{M},\mathcal{P},\mathcal{J}} F(U,\mathcal{M},\mathcal{P},\mathcal{J})\;\mathop = \limits^{\text{3}} \\ \mathop = \limits^{\text{3}} \;\mathop {min}\limits_U \left\{ {\mathop {min}\limits_{\mathcal{M},\mathcal{P},\mathcal{J}} F(U,\mathcal{M},\mathcal{P},\mathcal{J}\,{\text{|}}\,U)} \right\}\;\mathop = \limits^4 \;\mathop {min}\limits_U \mathop {min}\limits_{\mathcal{M},\mathcal{P},\mathcal{J}} G(\mathcal{M},\mathcal{P},\mathcal{J})\;\mathop = \limits^{\text{5}} \;\mathop {min}\limits_U {{G}_{U}}\;\mathop = \limits^{\text{6}} \;{{F}_{A}}. \\ \end{gathered} $

В этой цепочке равенство 1 – определение оптимального решения pадачи 1. Равенство 2 следует из вида целевой функции (1.1). Равенство 3 очевидно. Равенство 4 следует из формулы (4.1), равенство 5 – определение оптимального решения задачи 2 при каждом фиксированном значении $U$. Наконец, равенство 6 следует из формулы (5.1) в пошаговой записи алгоритма $\mathcal{A}$.

В соответствии с уверждением 2, время работы алгоритма ${{\mathcal{A}}_{1}}$ оценивается величиной $\mathcal{O}(T_{{{\text{max}}}}^{3}N)$. На шаге 1 алгоритма $\mathcal{A}$ требуется $K$ раз выполнить алгоритм ${{\mathcal{A}}_{1}}$, поэтому трудоемкость шага 1 алгоритма $\mathcal{A}$ есть величина $\mathcal{O}(KT_{{{\text{max}}}}^{3}N)$. Шаг 2, очевидно, требует $\mathcal{O}(K)$ операций, а шаг 3 – константного числа операций. Суммируя временные затраты на всех шагах алгоритма, получаем итоговую трудоемкость алгоритма $\mathcal{A}$. Теорема доказана.

Замечание 1. Если ${{T}_{{{\text{max}}}}}$ является частью входа задачи, а мощность $K$ алфавита – фиксированный параметр, то время работы алгоритма равно $\mathcal{O}({{N}^{4}})$. Таким образом, алгоритм $\mathcal{A}$ является полиномиальным.

6. ПРИМЕРЫ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ

Для иллюстрации работоспособности алгоритма приведем два примера обработки последовательностей (временных рядов), которые можно интерпретировать как квазипериодические последовательности флуктуирующих ECG-подобных (фиг. 4) и PPG-подобных (фиг. 5) импульсов при наличии аддитивной помехи. В действительности, с математической точки зрения не имеет значения, какие именно последовательности включены в состав алфавита $W$. Выбор именно ECG и PPG-подобных сигналов обусловлен желанием проиллюстрировать применимость алгоритма для биомедицинских приложений. Форма этих импульсов, а так же характерные участки и значимые точки описаны экспертами, см., например, [6]–[9]. На фиг. 4 (ECG-импульс) этими участками являются P-Q-R-S-T-U, на фиг. 5 (PPG-импульс) – систолическая точка, дикротическая точка, диастолическая точка. Раскраска фигур маркирует характерные участки или интервалы между значимыми точками.

Фиг. 4.

Пример 1. Распознавание ECG-подобной последовательности.

Фиг. 5.

Пример 2. Распознавание PPG-подобной последовательности.

Рассмотрим подробнее фиг. 4, иллюстрирующую помехоустойчивую обработку ECG-подобного сигнала. В верхнем ряду изображен алфавит. Ниже слева изображена последовательность $U$ – одна из последовательностей из алфавита. Справа от нее изображена программно-сгенерированная последовательность – квазипериодическая последовательность, порожденная флуктуирующими повторами импульса $U$. Под модельной последовательностью изображена последовательность $Y$ – результат сложения модельной последовательности и последовательности независимых одинаково распределенных гауссовских случайных величин с нулевым математическим ожиданием. Важно отметить, что только $Y$ и $W$ являются входными данными для алгоритма. Последовательность $U$ и модельная последовательность $X$, изображенные во втором ряду фигуры, приведены для иллюстрации; эти данные недоступны.

В нижней части фигуры представлены результаты работы алгоритма $\mathcal{A}$: слева последовательности ${{U}_{A}}$, которая является оценкой для модельной, ненаблюдаемой последовательности. Справа изображены компоненты последовательности ${{X}_{A}}$, восстановленной по формулам (1.7) и (1.8) с использованием четырех наборов ${{U}_{A}}$, ${{\mathcal{M}}_{A}}$, ${{\mathcal{P}}_{A}}$ и ${{\mathcal{J}}_{A}}$, полученных в результате работы алгоритма $\mathcal{A}$. Данные на фиг. 4 получены при $K = 3$, $q(U) = 203$, $U \in W$, ${{T}_{{{\text{max}}}}} = 370$, $N = 1800$, максимальная амплитуда сигнала – $128$, уровень шума $\sigma = 35$.

Фигура 5 показывает результаты обработки PPG-подобного сигнала. Эта фигура имеет такую же структуру, как и фиг. 4. Данные на фиг. 5 получены при $K = 3$, $q(U) = 120$, $U \in W$, ${{T}_{{{\text{max}}}}} = 240$, $N = 1000$, максимальная амплитуда сигнала – 128, уровень шума $\sigma = 35$.

Приведенные примеры показывают, что построенный алгоритм позволяет с вполне приемлемым качеством обрабатывать данные в виде квазипериодических последовательностей флуктуирующих импульсов. Во-первых, неизвестная последовательность $U$ и определенная алгоритмом последовательность ${{U}_{A}}$ совпадают в обоих примерах (распознавание произведено корректно). Во-вторых, визуальное сравнение графиков ненаблюдаемой последовательности $X$ и восстановленной последовательности ${{X}_{A}}$ демонстрирует лишь незначительные отклонения одного графика от другого и практически точное совпадение всех характерных точек.

ЗАКЛЮЧЕНИЕ

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

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

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

Список литературы

Kel’manov A.V., Khamidullin S.A., Okol’nishnikova L.V. Recognition of a quasiperiodic sequence containing identical subsequences-fragments // Pattern Recognition and Image Analysis. 2004. V. 14. № 1. P. 72–83.
Kel’manov A.V., Khamidullin S.A. Recognizing a quasiperiodic sequence composed of a given number of identical subsequences // Pattern Recognition and Image Analysis. 2000. V. 10. № 1. P. 127–142.
Kel’manov A.V., Khamidullin S.A., Mikhailova L.V., Ruzankin P.S. Polynomial-time solvability of one optimization problem induced by processing and analyzing quasiperiodic ECG and PPG signals // Communications in Computer and Information Science. 2020. V. CCIS 1145. P. 88–101.
Kel’manov A.V., Khamidullin S.A. Algorithm of recognition of a quasiperiodic sequence composed of a given number of truncated subsequences // Pattern Recognition and Image Analysis. 2001. V. 11. № 1. P. 43–46.
Кельманов А.В., Хамидуллин С.А. Распознавание числовой последовательности по фрагментам квазипериодически повторяющейся эталонной последовательности // Сиб. ж. индустр. матем. 2004. Т. 7. № 2. С. 68–87.
Rajni R., Kaur I. Electrocardiogram signal analysis – an overview // Int. J. Comput. Appl. 2013. V. 84. № 7. P. 22–25.
Al-Ani M.S. ECG Waveform Classification Based on P-QRS-T Wave Recognition // UHD Journal of Science and Technology. 2018. V. 2. № 2.
Shelley K., Shelley S. Pulse Oximeter Waveform: Photoelectric Plethysmography. In: Carol Lake, Hines R., Blitt C. (eds). Clinical Monitoring. W.B. Saunders Company. 2001. P. 420–428.
Elgendi M. On the analysis of fingertip photoplethysmogram signals // Current Cardiology Reviews. 2012. V. 8. № 1. P. 14–25.

Дополнительные материалы отсутствуют.

Инструменты

следующая статья выпуска предыдущая статья выпуска содержание выпуска

Журнал вычислительной математики и математической физики

Архивы выпусков Информация о журнале Отправить рукопись в журнал