Динамическая деформация времени для геологических временных рядов, Matlab

Question

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

здесь dtw в matlab:

function [Dist,D,k,w]=dtw()
%Dynamic Time Warping Algorithm
%Dist is unnormalized distance between t and r
%D is the accumulated distance matrix
%k is the normalizing factor
%w is the optimal path
%t is the vector you are testing against
%r is the vector you are testing
[t,r,x1,x2]=randomtestdata();
[rows,N]=size(t);
[rows,M]=size(r);
%for n=1:N
%    for m=1:M
%        d(n,m)=(t(n)-r(m))^2;
%    end
%end
d=(repmat(t(:),1,M)-repmat(r(:)',N,1)).^2; %this replaces the nested for loops from         above Thanks Georg Schmitz 

D=zeros(size(d));
D(1,1)=d(1,1);

for n=2:N
    D(n,1)=d(n,1)+D(n-1,1);
end
for m=2:M
    D(1,m)=d(1,m)+D(1,m-1);
end
for n=2:N
    for m=2:M
        D(n,m)=d(n,m)+min([D(n-1,m),D(n-1,m-1),D(n,m-1)]);
    end
end

Dist=D(N,M);
n=N;
m=M;
k=1;
w=[];
w(1,:)=[N,M];
while ((n+m)~=2)
    if (n-1)==0
        m=m-1;
    elseif (m-1)==0
        n=n-1;
    else 
      [values,number]=min([D(n-1,m),D(n,m-1),D(n-1,m-1)]);
      switch number
      case 1
        n=n-1;
      case 2
        m=m-1;
      case 3
        n=n-1;
        m=m-1;
      end
  end
    k=k+1;
    w=cat(1,w,[n,m]);
end
w=flipud(w)

%w is a matrix that looks like this:

%    1 1
%    1 2
%    2 2
%    3 3
%    3 4
%    3 5
%    4 5
%    5 6
%    6 6

, так что это говорит о том, что как первая, так и вторая точки второго вектора должны быть сопоставлены с первой точкой первого вектора,т.е. 1 1 1 2 и что пятая и шестая точки на первом векторе должны быть сопоставлены со вторым вектором в точке шесть.и т.д., поэтому w содержит координаты x деформированных данных.

Обычно я могу сказать

X1=w(:,1);
X2=w(:,2);
for i=1:numel(reference vector)
 Y1(i)=reference vector(X1(i));
 Y2(i)=test vector(X2(i));
end

, но мне не нужно растягивать опорный вектор, поэтому мне нужно использоватьповторяет в X1, чтобы знать, как сжать Y2, и повторяет в X2, чтобы знать, как растянуть Y2, вместо того, чтобы использовать повторы в X1 для растяжения Y1 и повторяет в X2, чтобы растянуть Y2.

Я пытался использовать метод findнайти повторы в X1 и X2, а затем усреднить (сжать) или интерполировать линейно (растянуть) по мере необходимости, но код стал очень сложным и трудным для отладки.

Было ли это действительно неясно?Мне было трудно объяснить эту проблему, но мне просто нужно знать, как взять w и создать Y2, который соответственно растягивается и сокращается.

Answer 1

Во-первых, вот DTW в Matlab, переведенный с псевдокода на wikipedia :

t = 0:.1:2*pi;
x0 = sin(t) + rand(size(t)) * .1;
x1 = sin(.9*t) + rand(size(t)) * .1;

figure
plot(t, x0, t, x1);
hold on

DTW = zeros(length(x0), length(x1));
DTW(1,:) = inf;
DTW(:,1) = inf;
DTW(1,1) = 0;

for i0 = 2:length(x0)
    for i1 = 2:length(x1)
        cost = abs(x0(i0) - x1(i1));
        DTW(i0, i1) = cost + min( [DTW(i0-1, i1) DTW(i0, i1-1) DTW(i0-1, i1-1)] );
    end
end

Если вы деформируете x_0 на x_1, x_1 на x_0 или деформируете их друг на друга, вы можете получить ответ из матрицы DTW. В вашем случае:

[cost, path] = min(DTW, [], 2);
plot(t, x1(path));
legend({'x_0', 'x_1', 'x_1 warped to x_0'});