用MATLAB实现DFA算法及报告

地理与生物信息学院2011/2012学年第二学期实验报告课程名称：医学信号处理实验名称：DFA算法实现班级学号：**********学生姓名：****指导教师：****日期：2012年6月2一、实验题目：Matlab软件处理生物医学数字信号二、实验要求：采用DFA算法实现三、实验内容：能使用Matlab软件编写、实现DFA算法并处理生物医学信号四、实验过程：（1）针对某个心电时间序列(总长度为N)，首先进行求和])([)(1kiaveBiBky其中，)(iB是第i个数据，而aveB是所分析的心电时间序列的平均值。这种求和步骤可以映射原始时间序列到一个自相似的过程。（2）接着，量度求和后的时间序列的垂直特征尺度。方法是：求和后的时间序列被分成等长为n的许多个小片段。画出每一个长度为n的小片段的最小二乘拟和直线（它用来代表该片段中的趋势）。直线段的y坐标标记为)(kyn。（3）其次，对求和后的时间序列进行去趋势，即在每一个片段中把)(ky减去局域趋势)(kyn。对一个给定的片段长度n，这种求和并且去趋势后的时间序列的涨落的特征尺寸可以由下式计算：NknkykyNnF12)]()([1)(（4）在所有时间尺度（片段大小）上重复上面的计算，可以得出)(nF和片段大小n的关系曲线。)(lognF对nlog的斜率决定尺度指数（自相似参量）。3五、主要实验程序：Clearall;A=textread('d:\nsrdb-16786.txt');N=1000;n=[48163264128];l=length(n);r=fix(rand(1)*10000);b(1:N)=A(r:(r-1+N),2)*100+2100;plot(1:N,b(1:N));c=b-mean(b);fori=1:Ny(i)=sum(c(1:i));endfork=1:l%k=1n=10k=2n=50...d=N/n(k);forj=1:dp=polAit((j-1)*n(k)+1:j*n(k),y((j-1)*n(k)+1:j*n(k)),1);y1(k,(j-1)*n(k)+1:j*n(k))=polyval(p,(j-1)*n(k)+1:j*n(k));endfigure,plot(1:N,y);holdon;plot(1:N,y1(k,1:N),'r');set(gca,'XGrid','on');%set(gca,'XTick',[0:n(k):N]);4title(['n=',num2str(n(k))]);f(k)=sqrt((sum((y(1:N)-y1(k,1:N)).^2))/N);endfigure,plot(log10(n(1:k)),log10(f(1:k)),'o');holdon;a=polAit(log10(n(1:k)),log10(f(1:k)),1);plot(log10(n(1:k)),polyval(a,log10(n(1:k))));title(['拟合直线斜率',num2str(a(1))]);原始数据曲线：n=4、8、16、32、64、12856拟合直线：7六、实验小结：DFA方法可以检测出一个非平稳时间序列的内在的自相似性，同时避免散杂地检测可明显看出的人为的外在趋势的自相似性。更适合应用于某些慢性趋势非平稳时间序列。通过此次实验，熟悉了DFA算法的Matlab实现，并运用其处理生物医学信号。