matlab__數(shù)據(jù)分析和統(tǒng)計(jì)資料

1、Good is good, but better carries it.精益求精，善益求善。matlab_數(shù)據(jù)分析和統(tǒng)計(jì)-數(shù)據(jù)分析和統(tǒng)計(jì)面向列的數(shù)據(jù)集這年頭似乎十分風(fēng)行”面向”這個(gè)詞，這兒故也套用，其英文為Column-OrientedDataSets,可理解為MatLab按列的存儲(chǔ)方式來(lái)分析數(shù)據(jù),下面是一個(gè)例子:TimeLocation1Location2Location301h001111902h007131103h0014172004h001113905h0043516906h0038467607h006113218608h007513518009h00388811510h0028365

2、511h0012121412h0018273013h0018192914h0017151815h0019364816h0032471017h0042659218h00576615119h0044559020h0011414525721h0035586822h0011121523h001391524h001097以上數(shù)據(jù)被保存在一個(gè)稱(chēng)為count.dat的文件中11119713111417201113943516938467661132186751351803888115283655121214182730181929171518193648324710426592576615144559011

3、4145257355868111215139151097下面，我們調(diào)入此文件，并看看文件的一些參數(shù)loadcount.datn,p=size(count)n=24p=3創(chuàng)建一個(gè)時(shí)間軸后，我們可以把圖畫(huà)出來(lái)：t=1:n;set(0,defaultaxeslinestyleorder,-|-|-.)set(0,defaultaxescolororder,000)plot(t,count),legend(Location1,Location2,Location3,0)xlabel(Time),ylabel(VehicleCount),gridon足以證明，以上是對(duì)個(gè)對(duì)象的次觀測(cè)基本數(shù)據(jù)分析函數(shù)（一定

4、注意是面向列的）繼續(xù)用上面的數(shù)據(jù)，其每列最大值均值及偏差分別為：mx=max(count)mu=mean(count)sigma=std(count)mx=114145257mu=32.000046.541765.5833sigma=25.370341.405768.0281重載函數(shù)，還可以定位出最大最小值的位置mx,indx=min(count)mx=797indx=22324試試看，你能看懂下面的命令是干什么的嗎？n,p=size(count)e=ones(n,1)x=counte*muHYPERLINKjavascript:void(null)點(diǎn)這看看答案！下面這句命令則找出了整個(gè)矩陣的

5、最小值:min(count(:)ans=7協(xié)方差及相關(guān)系數(shù)下面,我們來(lái)看看第一列的方差:cov(count(:,1)ans=643.6522cov()函數(shù)作用于矩陣,則會(huì)計(jì)算其協(xié)方差矩陣.corrcoef()用于計(jì)算相關(guān)系數(shù),如:corrcoef(count)ans=1.00000.93310.95990.93311.00000.95530.95990.95531.0000數(shù)據(jù)的預(yù)處理未知數(shù)據(jù)NaN(NotaNumber-不是一個(gè)數(shù))被定義為未經(jīng)定義的算式的結(jié)果,如0/0.在處理數(shù)據(jù)中,NaN常用來(lái)表示未知數(shù)據(jù)或未能獲得的數(shù)據(jù).所有與NaN有關(guān)的運(yùn)算其結(jié)果都是NaN.a=magic(3);a(

6、2,2)=NaNa=8163NaN7492sum(a)ans=15NaN15在做統(tǒng)計(jì)時(shí)，常需要將NaN轉(zhuǎn)化為可計(jì)算的數(shù)字或去掉,以下是幾種方法:注:判斷一個(gè)值是否為NaN,只能用isnan(),而不可用x=NaN;i=find(isnan(x);x=x(i)先找出值不是NaN的項(xiàng)的下標(biāo),將這些元素保留x=x(find(isnan(x)同上,去掉NaNx=x(isnan(x);更快的做法x(isnan(x)=;消掉NaNX(any(isnan(X),:)=;把含有NaN的行都去掉用此法可以從數(shù)據(jù)中去掉不相關(guān)的數(shù)據(jù),看看下面的命令是干什么用的:mu=mean(count);sigma=std(co

7、unt);n,p=size(count)outliers=abs(countmu(ones(n,1),:)3*sigma(ones(n,1),:);nout=sum(outliers)nout=100count(any(outliers),:)=;HYPERLINKjavascript:void(null)點(diǎn)這看看答案回歸與曲線(xiàn)擬合我們經(jīng)常需要把觀測(cè)到的數(shù)據(jù)表達(dá)為函數(shù)，假如有如下的對(duì)時(shí)間的觀測(cè)：t=1.62.3;y=0.50.851.40;plot(t,y,o),gridon多項(xiàng)式回歸由圖可以看出應(yīng)該可以用多項(xiàng)式來(lái)表達(dá)：y=a0+a1*t+a2*t2系數(shù)a0,a1,a2可以由最小平方擬合來(lái)確定

8、，這一步可由反除號(hào)來(lái)完成解下面的三元方程組可得：X=ones(size(t)tt.2X=1.0000001.00000.30000.09001.00000.80000.64001.00001.10001.21001.00001.60002.56001.00002.30005.2900a=Xya=0.53180.91910.2387a即為待求的系數(shù)，畫(huà)圖比較可得T=(0:0.1:2.5);Y=ones(size(T)TT.2*a;plot(T,Y,t,y,o,),gridon結(jié)果令人失望，但我們可以增加階數(shù)來(lái)提高精確度，但更明智的選擇是用別的方法線(xiàn)性參數(shù)回歸形如：y=a0+a1*exp(-t)+

9、a2*t*exp(-t)計(jì)算方法同上：X=ones(size(t)exp(t)t.*exp(t);a=Xya=1.39740.89880.4097T=(0:0.1:2.5);Y=ones(size(T)exp(T)T.exp(T)*a;plot(T,Y,t,y,o),gridon看起來(lái)是不是好多了！例子研究：曲線(xiàn)擬合下面我們以美國(guó)人口普查的數(shù)據(jù)來(lái)研究一下有關(guān)曲線(xiàn)擬合的問(wèn)題（MatLab是別人的,教學(xué)文檔是別人的,例子也是別人的，我只是一個(gè)翻譯而已）loadcensus這樣我們得到了兩個(gè)變量，cdate是1790至1990年的時(shí)間列向量(10年一次),pop是相應(yīng)人口數(shù)列向量.上一小節(jié)所講的多項(xiàng)

10、式擬合可以用函數(shù)polyfit()來(lái)完成,數(shù)字指明了階數(shù)p=polyfit(cdate,pop,4)Warning:Matrixisclosetosingularorbadlyscaled.Resultsmaybeinaccurate.RCOND=5.429790e20p=1.0e+05*0.00000.00000.00000.01266.0020產(chǎn)生警告的原因是計(jì)算中的cdata值太大，在計(jì)算中的Vandermonde行列式使變換產(chǎn)生了問(wèn)題，解決的方法之一是使數(shù)據(jù)標(biāo)準(zhǔn)化預(yù)處理：標(biāo)準(zhǔn)化數(shù)據(jù)數(shù)據(jù)的標(biāo)準(zhǔn)化是對(duì)數(shù)據(jù)進(jìn)行縮放，以使以后的計(jì)算能更加精確，一種方法是使之成為0均值：sdate=(cdate

11、mean(cdate)./std(cdate)現(xiàn)在再進(jìn)行曲線(xiàn)擬合就沒(méi)事了！p=polyfit(sdate,pop,4)p=0.70470.921023.470673.859862.2285pop4=polyval(p,sdate);plot(cdate,pop4,cdate,pop,+),gridon在上面的數(shù)據(jù)標(biāo)準(zhǔn)化中，也可以有別的算法，如令年的人口數(shù)為余量分析p1=polyfit(sdate,pop,1);pop1=polyval(p1,sdate);plot(cdate,pop1,cdate,pop,+)res1=poppop1;figure,plot(cdate,res1,+)p=po

12、lyfit(sdate,pop,2);pop2=polyval(p,sdate);plot(cdate,pop2,cdate,pop,+)res2=poppop2;figure,plot(cdate,res2,+)p=polyfit(sdate,pop,4);pop4=polyval(p,sdate);plot(cdate,pop4,cdate,pop,+)res4=poppop4;figure,plot(cdate,res4,+)可以看出，多項(xiàng)式擬合即使提高了階次也無(wú)法達(dá)到令人滿(mǎn)意的結(jié)果指數(shù)擬合從人口增長(zhǎng)圖可以發(fā)現(xiàn)人數(shù)的增長(zhǎng)基本是呈指數(shù)增加的，因此我們可以用年份的對(duì)數(shù)來(lái)進(jìn)行擬合，這兒，年數(shù)是

13、標(biāo)準(zhǔn)化后的！logp1=polyfit(sdate,log10(pop),1);logpred1=10.polyval(logp1,sdate);semilogy(cdate,logpred1,cdate,pop,+);gridonlogres2=log10(pop)polyval(logp2,sdate);plot(cdate,logres2,+)上面的圖不令人滿(mǎn)意，下面，我們用二階的對(duì)數(shù)分析：logp2=polyfit(sdate,log10(pop),2);logpred2=10.polyval(logp2,sdate);semilogy(cdate,logpred2,cdate,pop

14、,+);gridonr=pop10.(polyval(logp2,sdate);plot(cdate,r,+)這種余量分析比多項(xiàng)式擬合的余量分析圖案要隨機(jī)的多（沒(méi)有很強(qiáng)的規(guī)律性），可以預(yù)見(jiàn)，隨著人數(shù)的增加，余糧所反映的不確定度也在增加，但總的說(shuō)來(lái)，這種擬合方式要強(qiáng)好多！誤差邊界誤差邊界常用來(lái)反映你所用的擬合方式是否適用于數(shù)據(jù)，為得到誤差邊界，只需在polyfit()中傳遞第二個(gè)參數(shù)，并將其送入polyval().下面是一個(gè)二階多項(xiàng)式擬合模型,年份已被標(biāo)準(zhǔn)化,下面的代碼用了2,對(duì)應(yīng)于95%的可置信度:p2,S2=polyfit(sdate,pop,2);pop2,del2=polyval(p2,

15、sdate,S2);plot(cdate,pop,+,cdate,pop2,g,cdate,pop2+2*del2,r:,.cdate,pop22*del2,r:),gridon差分方程和濾波MatLab中的差分和濾波基本都是對(duì)向量而言的，向量則是存儲(chǔ)取樣信號(hào)或序列的函數(shù)y=filter(b,a,x)將用a,b描述的濾波器處理向量x,然后將其存儲(chǔ)在向量y中,filter()函數(shù)可看為是一差分方程a1y(n)=b1*x(1)+b2*x(2)+.-a2*y(2)-.如有以下差分方程:y(n)=1/4*x(n)+1/4*x(n-1)+1/4*x(n-2)+1/4*x(n-3),則a=1;b=1/41

16、/41/41/4;我們載入數(shù)據(jù),取其第一列,并計(jì)算有:loadcount.datx=count(:,1);y=filter(b,a,x);t=1:length(x);plot(t,x,.,t,y,),gridonlegend(OriginalData,SmoothedData,2)實(shí)現(xiàn)所表示的就是濾波后的數(shù)據(jù),它代表了4小時(shí)的平均車(chē)流量MatLab的信號(hào)處理工具箱中提供了很多用來(lái)濾波的函數(shù),可用來(lái)處理實(shí)際問(wèn)題!快速傅立葉變換(FFT)傅立葉變換能把信號(hào)按正弦展開(kāi)成不同的頻率值,對(duì)于取樣信號(hào),用的是離散傅立葉變換.FFT是離散傅立葉變換的一種高速算法,在信號(hào)和圖像處理中有極大的用處!fft離散傅

17、立葉變換fft2二維離散傅立葉變換fftnn維離散傅立葉變換ifft離散傅立葉反變換ifft2二維離散傅立葉反變換ifftnn維離散傅立葉反變換abs幅度angle相角unwrap相位按弧度展開(kāi),大于的變換為2的補(bǔ)角fftshift把零隊(duì)列移至功率譜中央cplxpair把數(shù)據(jù)排成復(fù)數(shù)對(duì)nextpow2下兩個(gè)更高的功率向量x的FFT可以這樣求:x=43791000y=fft(x)y=6.000011.48532.7574i2.000012.0000i5.4853+11.2426i18.00005.485311.2426i2.0000+12.0000i11.4853+2.7574ix雖然是實(shí)數(shù),但

18、y是復(fù)數(shù),其中,第一個(gè)是因?yàn)樗浅?shù)相加的結(jié)果,第五個(gè)則對(duì)應(yīng)于奈奎斯特頻率,后三個(gè)數(shù)是由于負(fù)頻率的影響,它們是前面三個(gè)數(shù)的共軛值!下面,讓我們來(lái)驗(yàn)證一下太陽(yáng)黑子活動(dòng)周期是11年!Wolfer數(shù)記錄了300年太陽(yáng)黑子的數(shù)量及大小:loadsunspot.datyear=sunspot(:,1);wolfer=sunspot(:,2);plot(year,wolfer)title(SunspotData)現(xiàn)在來(lái)看看其FFT:Y=fft(wolfer);Y的幅度是功率譜,畫(huà)出功率譜和頻率的對(duì)應(yīng)關(guān)系就得出了周期圖,去掉第一點(diǎn),因?yàn)樗皇撬袛?shù)據(jù)的和,畫(huà)圖有:N=length(Y);Y(1)=;power=abs(Y(1:N/2).2;nyquist=1/2;freq=(1:N/2)/(N/2)*nyquist;plot(freq,power),gridonxlabel(cycles/year)title(Periodogram)上面的圖看起來(lái)不大方便,下面我們畫(huà)出頻譜-周期圖period=1./freq;plot(period,power),axis(04002e7),gridonylabel(Power)xlab