最重點的事情在本篇的最后，你也可以直接看末尾——medfilt優(yōu)先級最高

數(shù)據(jù)預(yù)處理，顧名思義，預(yù)處理操作有很多，總的來說，它的一個形象的起搏就是“胃腸前端的咀嚼”，較為復(fù)雜的預(yù)處理和正式處理的界限會模糊——對于消化系統(tǒng)，胃其實也算腸的預(yù)處理前端，但它比咀嚼要復(fù)雜的多。

我們在咀嚼時，會把硌牙的異物吐出來已免其進(jìn)入腸胃——剔除異常值，為此需要判定硌牙的程度——置信

拉依達(dá)方法就是高中數(shù)學(xué)所學(xué)的3西格瑪原則，大于3西格瑪原則的我們認(rèn)為它不可信，也就是認(rèn)為它是假的，就把它剔除。話說代碼實例略是幾個意思？考題不會出現(xiàn)？或者說，這事干脆不用MATLAB，用Excel/WPS表格就行?——但我用后者發(fā)現(xiàn)不會弄……還是用MATLAB自己寫一個吧~

我們現(xiàn)在來把異常值踢了：

行列兩重循環(huán)遍歷所有元素，將異常值置0

這個wn叫肖偉勒系數(shù)，它隨著測量次數(shù)n的增加而緩慢增加，

Sx是數(shù)據(jù)xi的標(biāo)準(zhǔn)差，拉依達(dá)方法的3西格瑪原則相當(dāng)于讓肖偉勒系數(shù)wn=3。

x=load('error.dat');%導(dǎo)入數(shù)據(jù)

n=length(x);%求得測量次數(shù)n

subplot(2,1,1);

plot(x,'o');

title('原始數(shù)據(jù)')

axis([0,n+1,min(x)-1,max(x)+1]);%axis函數(shù)設(shè)置圖的橫縱視野范圍比數(shù)據(jù)范圍大1

w=1+0.4*log(n);%要求不嚴(yán)格時的近似公式，話說我不清楚怎么算是嚴(yán)格

yichang = abs(x-mean(x)) > w*std(x);

% 若用拉依達(dá)方法，把w改成3即可，但本組數(shù)據(jù)將不能成功剔除異常值。

x(yichang)=[];

save errornew.dat x -ASCII

subplot(2,1,2);

plot(x,'rs');

title('異常值剔除后數(shù)據(jù)');

axis([0,n+1,min(x)-1,max(x)+1]);

與拉依達(dá)和肖維勒方法相比，下邊這種方法風(fēng)格差異較大：一階差分法

看起來這個方法似乎用不太上的樣子……得實時采集，還得有較強(qiáng)的線性，講義上提了一嘴，然后就話鋒一轉(zhuǎn)：

對于平滑處理，形象的起搏是“給食物剝皮”——噪聲就好比皮殼，它并不與“腸胃”親和。

如圖，我們傾向于認(rèn)為一段數(shù)據(jù)是趨于穩(wěn)定的，而突出的毛刺被置信為需要濾去噪聲。相比于剔除異常數(shù)據(jù)，平滑處理不僅是剔除了毛刺，還填充了一個我們的期望數(shù)據(jù)：

先建立一個函數(shù)文件Y，

% 建立“n點單純移動平均”的濾波函數(shù)

% 注意函數(shù)要單獨保存為與函數(shù)名同名的.m文件

function Y=smooth_data(y,n)??

m=length(y);

j=1;

for i=(n-1)/2+1:(m-(n-1)/2)

????p=i-(n-1)/2;

????q=i+(n-1)/2;

????Y(j)=sum(y(p:q))/n;

????j=j+1;

end

end

再在命令行窗口輸入：

% 主程序

clc

clear

t=-15:0.5:15;

n=length(t);

Y=5./(1+t.^2);???% 原始測試數(shù)據(jù)

y=Y+(0.5-rand(1,n));??% 給測試數(shù)據(jù)加上噪聲干擾

y1=smooth_data(y,9);??% 調(diào)用函數(shù)作9點濾波處理

plot(1:n,Y,1:n,y,'-o',5:n-4,y1,'-*');

legend('無噪聲','含噪聲','9點平滑后');

這個方法對于2n+1的各個點是平權(quán)的，當(dāng)然按鄰域影響觀，越近的應(yīng)該有越高的權(quán)重，于是我們還想不平權(quán)——加權(quán)。

權(quán)重系數(shù)可以采用最小二乘原理，使平滑后的數(shù)據(jù)以最小均方差逼近原始數(shù)據(jù)。即令：

話說，這操作求了個導(dǎo)，比賽時能用得上么……也沒給個代碼示例……或者這些原理干脆就是屁話，直接調(diào)用庫函數(shù)拉倒了~

比賽優(yōu)先級最高的應(yīng)該是smooths函數(shù)：

因為它給了個示例是讀取240條股市數(shù)據(jù)，并進(jìn)行平滑處理：

以下內(nèi)容copy講義：

p0=x(1:240,1)'; % 用開盤價所在列的前240條數(shù)據(jù)，注意若不轉(zhuǎn)置可能導(dǎo)致后面處理結(jié)果異常

subplot(2,2,1);

plot(p0,'k','LineWidth',1.5);??% 繪制平滑后曲線圖，黑色實線，線寬1.5

xlabel('觀測序號');

ylabel('股市日開盤價');

axis([0 250 1000 1400]);

p1 = smoothts(p0,'b',30); ??% 用盒子法平滑數(shù)據(jù)，窗寬為30

subplot(2,2,2);

plot(p0,'.'); ?% ?繪制日開盤價散點圖

plot(p0,'.','markersize',3); 可以改變點的大小

hold on

plot(p1,'k','LineWidth',1.5); ?

xlabel('觀測序號');

ylabel('盒子法');

legend('原始散點','平滑曲線','location','northwest');

axis([0 250 1000 1400]);

p2 = smoothts(p0,'g',30);% 高斯窗方法，窗寬為30，標(biāo)準(zhǔn)差為默認(rèn)值0.65

subplot(2,2,3);

plot(p0,'.');

hold on

plot(p2,'k','LineWidth',1.5); ??

xlabel('觀測序號'); ylabel('高斯窗方法');

legend('原始散點','平滑曲線','location','northwest');

axis([0 250 1000 1400]);

p3 = smoothts(p0,'e',30);??% 用指數(shù)法平滑數(shù)據(jù)，窗寬為30

subplot(2,2,4);

plot(p0,'.'); ?

hold on

plot(p3,'k','LineWidth',1.5);

xlabel('觀測序號'); ylabel('指數(shù)方法');

legend('原始散點','平滑曲線','location','northwest');

axis([0 250 1000 1400]);grid on

title('9點rloess平滑');

運行結(jié)果：

我個人感覺高斯窗是最為我們常用的，如圖它的平滑處理最貼近原始數(shù)據(jù)，這讓我想起了超分辨率熒光顯微技術(shù)中的STORM，用高斯擬合重建物像。

5、?用medfilt1函數(shù)（一維中值濾波）

medfilt顧名思義，med（medium）表示“中”，filt（filter）表示“過濾”，該篇文檔首先介紹的是“均值濾波”的平滑處理，用局部平均數(shù)代替噪聲值，現(xiàn)在的中值濾波就是用局部中位數(shù)代替噪聲值。

調(diào)用格式：

y = medfilt1(x,n)%由此，該函數(shù)的兩個參量是輸入信號和選定的取中值的窗口寬度

y = medfilt1(x,n,blksz)

y = medfilt1(x,n,blksz,dim)

例：?產(chǎn)生一列正弦波信號，加入噪聲信號，然后調(diào)用medfilt1函數(shù)對加入噪聲的正弦波進(jìn)行濾波（平滑處理）。

代碼：

t = linspace(0,4*pi,500)';? ?% 產(chǎn)生一個從0到4*pi的向量，長度為500個點，當(dāng)t軸用

y = 100*sin(t); ??% 產(chǎn)生正弦波信號

noise = normrnd(0,15,500,1);??% 產(chǎn)生500行1列的服從N(0,15)分布的隨機(jī)數(shù)，作為噪聲信號

%%話說講義上也不知咋回事，寫了個“服從N(0,152)分布”，盲猜是九鍵鍵盤5和2相鄰誤打了

y = y + noise; ??% 將正弦波信號疊加一個噪聲信號

subplot(2,1,1);

plot(t,y); ???

xlabel('時間');??ylabel('加噪聲的正弦波');

yy = medfilt1(y,30); ?%?調(diào)用medfilt1對加噪正弦波信號y進(jìn)行中值濾波、繪制窗寬為30的波形圖%講義上這個地方句子冗余，我給它改成這樣了

subplot(2,1,2);

plot(t,y,'b:'); ??% b:表示藍(lán)色虛線 ?

hold on

plot(t,yy,'k','LineWidth',2); % 繪制平滑后曲線，黑色實線，線寬2

xlabel('時間');??ylabel('中值濾波');??legend('加噪波形','平滑后波形');

運行結(jié)果：

yy還是500個點，沒有缺損

yy的每一個點都是局部取中位數(shù)，這個平滑作用是很強(qiáng)大的。

這個圖看起來效果上比smooth更好，這是因為——

“中值濾波是圖像處理中的一個常用步驟，它對于斑點噪聲（en:speckle noise）和椒鹽噪聲（en:salt-and-pepper noise）來說尤其有用。保存邊緣的特性使它在不希望出現(xiàn)邊緣模糊的場合也很有用。當(dāng)要求在降低噪聲的同時要求保持邊緣，中值濾波較卷積有更好的效果”——引用自https://blog.csdn.net/u010849228/article/details/48216919

我覺得，medfilt的優(yōu)先級最高