Есть два ряда чисел X(i), Y(i), i = { 0,...,n }
Требуется определить (хотя бы приблизительно) меру близости этих рядов с точностью до масштабного коэффициента их индексов, т.е. такую величину
Z = (Min по лямбда) ( (Сумма по i)( X(i)-Y(i*лямбда) )^2 )
Иными словами, вопрос — насколько близко можно свести ряды один к другому, растягивая или сжимая шаг по i у одного из них. Лямбда не обязательно целая, при этом Y(i*лямбда) понимается как величина, интерполированая по заданным для целых i значениям ряда.
Задачка взята не из теории. Поможите плиз любым советом. Бьюсь уже 2 недели.

Здравствуйте, Egen, Вы писали:

E>Есть два ряда чисел X(i), Y(i), i = { 0,...,n }
E>Требуется определить (хотя бы приблизительно) меру близости этих рядов с точностью до масштабного коэффициента их индексов, т.е. такую величину
E>Z = (Min по лямбда) ( (Сумма по i)( X(i)-Y(i*лямбда) )^2 )

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

Здравствуйте, Egen, Вы писали:

E>Есть два ряда чисел X(i), Y(i), i = { 0,...,n }
E>Требуется определить
E>Z = (Min по лямбда) ( (Сумма по i)( X(i)-Y(i*лямбда) )^2 )

Поскольку Y(i*лямбда) вычисляется, то, следовательно, имеется функция
z(Lambda), минимум которой (по Lambda) дает решение.
Смотрим программы минимизации функций одной переменной

http://www.srcc.msu.su/num_anal/lib_na/cat/cat122.htm

и выбираем по вкусу.
Константин.

Здравствуйте, Egen, Вы писали:

E>Есть два ряда чисел X(i), Y(i), i = { 0,...,n }
E>Требуется определить (хотя бы приблизительно) меру близости этих рядов с точностью до масштабного коэффициента их индексов, т.е. такую величину
E>Z = (Min по лямбда) ( (Сумма по i)( X(i)-Y(i*лямбда) )^2 )
E>Иными словами, вопрос — насколько близко можно свести ряды один к другому, растягивая или сжимая шаг по i у одного из них. Лямбда не обязательно целая, при этом Y(i*лямбда) понимается как величина, интерполированая по заданным для целых i значениям ряда.

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

            +inf
 F(T,f) = Integral X(t) * Y(t + T) * exp(-j*2*pi*t*F) * dt
            -inf

при фиксированном T = 0. (Извиняюсь, давно это было, поэтому мог и напутать )

Эффективно вычисляется с помощью быстрого преобразования Фурье (БПФ).

Здравствуйте, mrhru, Вы писали:

M>Если ряд Y(i) можно интерпретировать как частотный сдвиг ряда X(i) (плюс шумы),
M>то это задача вычисления свёртки в частотной области, т.е. вычисление максимума функции неопределённости
M>

M>            +inf
M> F(T,f) = Integral X(t) * Y(t + T) * exp(-j*2*pi*t*F) * dt
M>            -inf
M>

M>при фиксированном T = 0.

Извиняюсь, может я тупею, но вроде с переходом в частотную область масштабирование не сводится к сдвигу. Пусть например X(i)=синусоиды(2Гц+4Гц), а Y(i)=синусоиды(4Гц+8Гц).
X и Y можно совместить масштабированием индексов:Y(i/2)=X(i), однако в частотной области не удается совместить пики {2,4} с пиками {4,8} простым сдвигом. Там задача остается задачей масштабирования, как и во временной области
Вот если бы шкала в частотной области была бы нелинейной (не 1,2,3,4..., а скажем 1,2,4,8,...), тогда другое дело...

Здравствуйте, Egen, Вы писали:

E>Здравствуйте, mrhru, Вы писали:

M>>Если ряд Y(i) можно интерпретировать как частотный сдвиг ряда X(i) (плюс шумы),
M>>то это задача вычисления свёртки в частотной области, т.е. вычисление максимума функции неопределённости
M>>

M>>            +inf
M>> F(T,f) = Integral X(t) * Y(t + T) * exp(-j*2*pi*t*F) * dt
M>>            -inf
M>>

M>>при фиксированном T = 0.

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

Сводится.

E>Пусть например X(i)=синусоиды(2Гц+4Гц), а Y(i)=синусоиды(4Гц+8Гц).
E>X и Y можно совместить масштабированием индексов:Y(i/2)=X(i), однако в частотной области не удается совместить пики {2,4} с пиками {4,8} простым сдвигом.
Там задача остается задачей масштабирования, как и во временной области

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

E>Вот если бы шкала в частотной области была бы нелинейной (не 1,2,3,4..., а скажем 1,2,4,8,...), тогда другое дело...

Имхо, нет. Маштабирование по времени Y(i/2)=X(i) как раз и означает линейный сдвиг по частоте.

Здравствуйте, mrhru, Вы писали:

...

M>Имхо, нет. Маштабирование по времени Y(i/2)=X(i) как раз и означает линейный сдвиг по частоте.

Извиняюсь, по-моему — всё наврал.

Здравствуйте, klovetski, Вы писали:

K>Поскольку Y(i*лямбда) вычисляется, то, следовательно, имеется функция
K>z(Lambda), минимум которой (по Lambda) дает решение.
K>Смотрим программы минимизации функций одной переменной
K>и выбираем по вкусу.
K>Константин.

можно конечно и численными методами. Только вот с большим количеством участков немонотонности X и Y функция z приобретает катастрофически большое число локальных минимумов. Решение потребует большого времени и не факт, что сойдется куда надо.
Я думал идти другим путем. Локально представлять ряды полиномами 2-й степени. Тогда задача сводится к нахождению Lambda, минимизирующей
Z=(Сумма по t)[(a*t*t + b*t + c) — (A*t*t*Lambda*Lambda + B*t*Lambda + C)]^2,
а это сводится к кубическому уравнению по Lambda. Без численных методов.
Одно плохо. Есть подозрение, что a,b,c и A,B,C (коэф.полиномов) следует определять по выборкам разной длины {X(i),i=0...n} и {Y(j),j=0...N}, где N/n=Lambda. Если этого не делать, полиномы будут отражать локальные свойства рядов для разных (несопоставимых) участков

Здравствуйте, Egen, Вы писали:

E>Есть два ряда чисел X(i), Y(i), i = { 0,...,n }
E>Требуется определить (хотя бы приблизительно) меру близости этих рядов с точностью до масштабного коэффициента их индексов, т.е. такую величину
E>Z = (Min по лямбда) ( (Сумма по i)( X(i)-Y(i*лямбда) )^2 )
E>Иными словами, вопрос — насколько близко можно свести ряды один к другому, растягивая или сжимая шаг по i у одного из них. Лямбда не обязательно целая, при этом Y(i*лямбда) понимается как величина, интерполированая по заданным для целых i значениям ряда.
E>Задачка взята не из теории. Поможите плиз любым советом. Бьюсь уже 2 недели.

Можно попробовать произвести замену переменной j=ln(i). Таким образом мы сведем масштабирование к сдвигу. Затем вычислить корреляцинную функцию, как предлагал mrhru (правда, имхо, там формула несколько неправильная). Поскольку шаг после замены получается неравномерный, то придется прибегнуть к интерполяции.

Проблему локальных максимумов это до конца не решит, но во всяком случае получим более-менее внятную фунуцию.

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

Здравствуйте, <Аноним>, Вы писали:

А>Спасибо. Но этот и другие численные методы не подойдут по критерию времени и надежности. Основная идея — пожертвовать точностью, но выиграть время.

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

Здравствуйте, Кодт, Вы писали:

К>А если сделать спектральный анализ (даже самый дубовый, найти первые 2-3 гармоники) и сопоставить их частоты у обоих рядов?

Подумаю. Правда много всяких "но"...:
1. Сопоставить частоты — а как их сопоставишь, если спектр ряда, растянутого относительно другого ряда тоже растягивается (а не сдвигается!).
2. Ряд не обязательно гармонический (начало и конец выборки могут не совпадать),
3. анализ по гармоникам — но каким? Надо ведь по каким-то главным...
Может действительно, выбрать частоты с наибольшей амплитудой для X и Y, их и сопоставлять?
Спасибо за идею.

Здравствуйте, MichaelP, Вы писали:

MP>Можно попробовать произвести замену переменной j=ln(i). Таким образом мы сведем масштабирование к сдвигу. Затем вычислить корреляцинную функцию, как предлагал mrhru (правда, имхо, там формула несколько неправильная). Поскольку шаг после замены получается неравномерный, то придется прибегнуть к интерполяции.

MP>Проблему локальных максимумов это до конца не решит, но во всяком случае получим более-менее внятную фунуцию.


Конечно же! Как это я сам не догадался!?
Большое спасибо.

Здравствуйте, Egen, Вы писали:

E>1. Сопоставить частоты — а как их сопоставишь, если спектр ряда, растянутого относительно другого ряда тоже растягивается (а не сдвигается!).

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

E>2. Ряд не обязательно гармонический (начало и конец выборки могут не совпадать),

Можно не Фурье, а Лапласа... Я, правда, не скажу навскидку, как сделать быстрое преобразование Лапласа.

Если последовательности "пляшут" вокруг некоторых прямых ax+b, то вообще делать нечего Сравниваем коэффициенты при x.

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

Здравствуйте, Egen, Вы писали:

E>
E>Конечно же! Как это я сам не догадался!?
E>Большое спасибо.

Забыл упомянуть об одном достоинстве этого метода. Если считать корреляционную ф-цию через FFT, ее легко можно отмасштабировать так, чтобы ее значение давало бы коэф-т корреляции для данного "сдвига".

К>Если последовательности "пляшут" вокруг некоторых прямых ax+b, то вообще делать нечего Сравниваем коэффициенты при x.

К>Постоянную составляющую (среднее арифметическое) исключаем.
К>Дальше, если есть некоторый коэффициент роста (та самая ax), то все очень просто становится (чтобы найти, просто считаем сумму ряда, делим на длину).
К>Если нет — находим спектр.

Даже больше скажу: Если последовательности "пляшут" вокруг полиномов axx+bx+c, то тоже почти делать нечего. Метод наименьших квадратов дает аналитическую запись для коэффициента масштабирования индекса через решение уравнения Кардано (третьей степени)
Lambda^3 + A*Lambda + B = 0, где A,B,C выражаются через коэффициенты a,b,c аппроксимирующих полиномов для X и Y .
Но с попыткой перейти к полиномам более высоких порядков ( 3,4,5...) при аппроксимации X и Y такой подход приводит к уравнениям 5,7,9 ... степени, для которых современная математика пока слаба . Пользоваться же полиномом 2-й степени для представления рядов X и Y зачастую недостаточно.

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

Здравствуйте, Egen, Вы писали:

...

E>В форум я обратился в надежде прояснить свой замыленный взгляд. Чувствовал, что существует более тонкий метод. Просто не приходило в голову перейти к логарифмической шкале в частотной области . Здесь мне об этом тонко намекнули, за что большое спасибо. Буду копать здесь

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

MP>Может я невнятно выразился , но не надо переходить к логарифмической шкале в частотной области. Зачем делать дополнительную работу? Надо просто перейти к логарифмической шкале по x и считать корреляцию. Другое дело, что корреляцию может оказаться проще считать через FFT.

Пардон.
Вначале я понял Вас именно так, потом зачем-то переработал совет в частотную область (виноват, запутался в формализации логарифмического преобразования для нулевой точки последовательностей). Вы правы. Ln — он и во временной области Ln.

Постоянно забываю входить под своим именем...

Здравствуйте, Аноним, Вы писали:

А>Пардон.
А>Вначале я понял Вас именно так, потом зачем-то переработал совет в частотную область (виноват, запутался в формализации логарифмического преобразования для нулевой точки последовательностей). Вы правы. Ln — он и во временной области Ln.

Я вечером эксперимент в Excel провел.

Результат:
Для функции cos(x) на шестнадцати точках точность определения лямбда (2.5) — 1% . При этом коэф-т корреляции (при этом значении лямбда) прктически 1.

Правда были в эксперименте были некоторые особенности:
1. Чтобы не мучиться с мнимой частью в FFT эмулировал DCT (искусственно делал ф-цию симметричной)
2. Т.к. не знал, как повлияет цикличность свертки, чтобы избежать цикличности "добивал" функции нулями.
3. Оценивал положение максимума корреляционной функции квадратичной интерполяцией.
4. Самое существенное. Для простоты не интерополировал функции, а вычислял их в нужных точках.