珠海网站建设哪家好,龙华做网站,惠州百度seo在哪,黄埔网站推广本文只是简单的介绍DTW算法的目的和实现。具体的DTW可以参考一下文献#xff1a;
离散序列的一致性度量方法#xff1a;动态时间规整#xff08;DTW#xff09; http://blog.csdn.net/liyuefeilong/article/details/45748399
动态时间归整/规整/弯曲(Dynamic time warpi…本文只是简单的介绍DTW算法的目的和实现。具体的DTW可以参考一下文献
离散序列的一致性度量方法动态时间规整DTW http://blog.csdn.net/liyuefeilong/article/details/45748399
动态时间归整/规整/弯曲(Dynamic time warping,DTW) http://www.cnblogs.com/flypiggy/p/3603192.html
【更新日志】2017-7-17
评论区一位同学指出文章有误笔者能力有限暂未推导出错误之处这里提供一个一般不会出错的理论介绍网址维基百科毕竟能上维基百科的都是大家公认的。DTW相关维基百科戳这里有兴趣同学可以在评论区一起讨论文章错误之处谢谢
DTW是干什么的 动态时间规整算法故名思议就是把两个代表同一个类型的事物的不同长度序列进行时间上的“对齐”。比如DTW最常用的地方语音识别中同一个字母由不同人发音长短肯定不一样把声音记录下来以后它的信号肯定是很相似的只是在时间上不太对整齐而已。所以我们需要用一个函数拉长或者缩短其中一个信号使得它们之间的误差达到最小。 再来看看运动捕捉比如当前有一个很快的拳击动作另外有一个未加标签的动作我想判断它是不是拳击动作那么就要计算这个未加标签的动作和已知的拳击动作的相似度。但是呢他们两个的动作长度不一样比如一个是100帧一个是200帧那么这样直接对每一帧进行对比计算到后面误差肯定很大那么我们把已知拳击动作的每一帧找到无标签的动作的对应帧中使得它们的距离最短。这样便可以计算出两个运动的相似度然后设定一个阈值满足阈值范围就把未知动作加上“拳击”标签。
DTW怎么计算
下面我们来总结一下DTW动态时间规整算法的简单的步骤
1. 首先肯定是已知两个或者多个序列但是都是两个两个的比较所以我们假设有两个序列A{a1,a2,a3,...,am} B{b1,b2,b3,....,bn}维度mn
2. 然后用欧式距离计算出每序列的每两点之间的距离D(ai,bj) 其中1≤i≤m1≤j≤n 画出下表 3. 接下来就是根据上图将最短路径找出来。从D(a1,a2)沿着某条路径到达D(am,bn)。找路径满足假如当前节点是D(ai,bj)那么下一个节点必须是在D(i1,j)D(i,j1)D(i1,j1)之间选择并且路径必须是最短的。计算的时候是按照动态规划的思想计算也就是说在计算到达第(i,j)个节点的最短路径时候考虑的是左上角也即第(i-1,j)、(i-1,j-1)、(i,j-1)这三个点到(i,j)的最短距离。
4. 接下来从最终的最短距离往回找到那条最佳的输出路径 从D(a1,b1)到D(am,bn)。他们的总和就是就是所需要的DTW距离
【注】如果不回溯路径直接在第3步的时候将左上角三个节点到下一个节点最短的点作为最优路径节点就是贪婪算法了。DTW是先计算起点到终点的最小值然后从这个最小值回溯回去看看这个最小值都经过了哪些节点。
举个栗子
已知两个列向量a[8 9 1]b[ 2 5 4 6]其中代表转置也就是把行向量转换为列向量了
求两个向量利用动态时间规整以后的最短距离
第一步计算对应点的欧式距离矩阵d【注意开根号】 6 3 4 27 4 5 31 4 3 5
第二步计算累加距离D【从6出发到达5的累加距离】 6 9 13 1513 10 14 1614 14 13 18
计算方法如下
D(1,1)d(1,1)6
D(1,2)D(1,1)d(1,2)9
...
D(2,2)min(D(1,2),D(1,1),D(2,1))d(2,2)6410
...
D(m,n)min(D(m-1,n),D(m-1,n-1),D(m,n-1))d(m,n)
网上有很多代码但是大部分代码有误比如网上流传比较多的错误版本http://www.cnblogs.com/luxiaoxun/archive/2013/05/09/3069036.html
正确的代码可以自己写挺简单但是我找了一个可视化的代码大家可以参考一下
dtw.m function [Dist,D,k,w,rw,tw]dtw(r,t,pflag)
%
% [Dist,D,k,w,rw,tw]dtw(r,t,pflag)
%
% 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
% rw is the warped r vector
% tw is the warped t vector
% pflag plot flag: 1 (yes), 0(no)
%
% Version comments:
% rw, tw and pflag added by Pau Mic[row,M]size(r); if (row M) Mrow; rr; end;
[row,N]size(t); if (row N) Nrow; tt; end;
dsqrt((repmat(r,1,N)-repmat(t,M,1)).^2); %this makes clear the above instruction Thanks Pau Mic
d
Dzeros(size(d));
D(1,1)d(1,1);for m2:MD(m,1)d(m,1)D(m-1,1);
end
for n2:ND(1,n)d(1,n)D(1,n-1);
end
for m2:Mfor n2:ND(m,n)d(m,n)min(D(m-1,n),min(D(m-1,n-1),D(m,n-1))); % this double MIn construction improves in 10-fold the Speed-up. Thanks Sven Mensingend
endDistD(M,N);
nN;
mM;
k1;
w[M N];
while ((nm)~2)if (n-1)0mm-1;elseif (m-1)0nn-1;else [values,number]min([D(m-1,n),D(m,n-1),D(m-1,n-1)]);switch numbercase 1mm-1;case 2nn-1;case 3mm-1;nn-1;endendkk1;w[m n; w]; % this replace the above sentence. Thanks Pau Mic
end% warped waves
rwr(w(:,1));
twt(w(:,2));%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if pflag% --- Accumulated distance matrix and optimal pathfigure(Name,DTW - Accumulated distance matrix and optimal path, NumberTitle,off);main1subplot(position,[0.19 0.19 0.67 0.79]);image(D);cmap contrast(D);colormap(cmap); % copper bone, gray imagesc(D);hold on;xw(:,1); yw(:,2);indfind(x1); x(ind)10.2;indfind(xM); x(ind)M-0.2;indfind(y1); y(ind)10.2;indfind(yN); y(ind)N-0.2;plot(y,x,-w, LineWidth,1);hold off;axis([1 N 1 M]);set(main1, FontSize,7, XTickLabel,, YTickLabel,);colorb1subplot(position,[0.88 0.19 0.05 0.79]);nticks8;ticksfloor(1:(size(cmap,1)-1)/(nticks-1):size(cmap,1));mxmax(max(D));mnmin(min(D));ticklabelsfloor(mn:(mx-mn)/(nticks-1):mx);colorbar(colorb1);set(colorb1, FontSize,7, YTick,ticks, YTickLabel,ticklabels);set(get(colorb1,YLabel), String,Distance, Rotation,-90, FontSize,7, VerticalAlignment,bottom);left1subplot(position,[0.07 0.19 0.10 0.79]);plot(r,M:-1:1,-b);set(left1, YTick,mod(M,10):10:M, YTickLabel,10*rem(M,10):-10:0)axis([min(r) 1.1*max(r) 1 M]);set(left1, FontSize,7);set(get(left1,YLabel), String,Samples, FontSize,7, Rotation,-90, VerticalAlignment,cap);set(get(left1,XLabel), String,Amp, FontSize,6, VerticalAlignment,cap);bottom1subplot(position,[0.19 0.07 0.67 0.10]);plot(t,-r);axis([1 N min(t) 1.1*max(t)]);set(bottom1, FontSize,7, YAxisLocation,right);set(get(bottom1,XLabel), String,Samples, FontSize,7, VerticalAlignment,middle);set(get(bottom1,YLabel), String,Amp, Rotation,-90, FontSize,6, VerticalAlignment,bottom);% --- Warped signalsfigure(Name,DTW - warped signals, NumberTitle,off);subplot(1,2,1);set(gca, FontSize,7);hold on;plot(r,-bx);plot(t,:r.);hold off;axis([1 max(M,N) min(min(r),min(t)) 1.1*max(max(r),max(t))]);grid;legend(signal 1,signal 2);title(Original signals);xlabel(Samples);ylabel(Amplitude);subplot(1,2,2);set(gca, FontSize,7);hold on;plot(rw,-bx);plot(tw,:r.);hold off;axis([1 k min(min([rw; tw])) 1.1*max(max([rw; tw]))]);grid;legend(signal 1,signal 2);title(Warped signals);xlabel(Samples);ylabel(Amplitude);end 测试代码 test.m clear
clc
a[8 9 1 9 6 1 3 5];
b[2 5 4 6 7 8 3 7 7 2];
[Dist,D,k,w,rw,tw] DTW(a,b,1);
fprintf(最短距离为%d\n,Dist)
fprintf(最优路径为)
w 测试结果
d 6 3 4 2 1 0 5 1 1 67 4 5 3 2 1 6 2 2 71 4 3 5 6 7 2 6 6 17 4 5 3 2 1 6 2 2 74 1 2 0 1 2 3 1 1 41 4 3 5 6 7 2 6 6 11 2 1 3 4 5 0 4 4 13 0 1 1 2 3 2 2 2 3最短距离为27
最优路径为
w 1 11 21 31 41 51 62 63 74 85 96 107 108 10 规整以后可视化如下 【更新日志2019-9-11】
为了验证结果正确性在python中也找到了一个工具库dtw
安装方法
pip install dtw
接收参数与返回参数
def dtw(x, y, dist, warp1, winf, s1.0):Computes Dynamic Time Warping (DTW) of two sequences.:param array x: N1*M array:param array y: N2*M array:param func dist: distance used as cost measure:param int warp: how many shifts are computed.:param int w: window size limiting the maximal distance between indices of matched entries |i,j|.:param float s: weight applied on off-diagonal moves of the path. As s gets larger, the warping path is increasingly biased towards the diagonalReturns the minimum distance, the cost matrix, the accumulated cost matrix, and the wrap path.
输入参数前面x和y是输入的两个序列点dist是计算距离的方法因为有很多距离计算的方法只不过通常使用欧式距离直接使用numpy里面的linalg.norm即可。
返回值这个d有点蒙但是如果想得到上面的27那个结果可以看accumulated cost matrix这个返回值wrap path就是两条路径的对应点。
调用方法(官方案例):
import numpy as np# We define two sequences x, y as numpy array
# where y is actually a sub-sequence from x
x np.array([2, 0, 1, 1, 2, 4, 2, 1, 2, 0]).reshape(-1, 1)
y np.array([1, 1, 2, 4, 2, 1, 2, 0]).reshape(-1, 1)from dtw import dtweuclidean_norm lambda x, y: np.abs(x - y)d, cost_matrix, acc_cost_matrix, path dtw(x, y, disteuclidean_norm)print(d)
本文代码和文章已经同步到微信公众号中公众号与本博客将持续同步更新运动捕捉、机器学习、深度学习、计算机视觉算法敬请关注
