MATLAB應用——資料分析與統計

數學建模是用數學方法解決各種實際問題的橋樑，它已經滲透到各個領域，而且發揮出越來越重要的作用。面對自然科學和工程應用中的難題，大部分人無從入手，而個別人卻能短時間內給出切實可行的解決方案，其差別往往在於

駕馭數學知識的能力不同

。現代計算機技術的應用不僅減少了計算錯誤，而且加強了數學應用者解決問題的能力。MATLAB是一款常用的資料處理軟體，為了更好的應用MATLAB軟體，我將整理好的MATLAB函式分享到今日頭條上，以利己利人查閱。

MATLAB提供的很多資料分析與統計函式都是面向列的，即矩陣中的每一列代表一個變數的多個觀測值，其列數對應於變數數，行數對應於測量點數。

max和min函式可求出資料的最大值和最小值，mean和std函式可求出資料的均值和標準差，sum和prod函式可求出資料元素和與資料元素積。例如，對MATLAB內含的某城市24小時的車流量資料count。dat可作分析：

load count。dat

mx=max（count）

mx = 114 145 257

mu=mean（count）

mu = 32。0000 46。5417 65。5833

sigma=std（count）

sigma = 25。3703 41。4057 68。0281

對有些函式還可給出位置，例如，在求出最小值的同時，可得到最小值所在的位置（行號）：

［mx，indx］=min（count）

mx = 7 9 7

indx = 2 23 24

1、協方差和相關係數

cov函式可以求出單個變數的協方差，而corrcoef函式可求出兩個變數之間的相關係數，例如：

cv=cov（count）

cv = 1。0e+003 *

0。6437 0。9802 1。6567

0。9802 1。7144 2。6908

1。6567 2。6908 4。6278

cr=corrcoef（count）

cr =

1。0000 0。9331 0。9599

0。9331 1。0000 0。9553

0。9599 0。9553 1。0000

2、資料預處理

在MATLAB中遇到超出範圍的資料時均用NaN （非數值） 表示，而且在任何運算中，只要包含NaN，就將它傳遞到結果中，因此在對資料進行分析前，應對資料中出現的NaN作剔除處理。例如：

a=［1 2 3；5 NaN 8；7 4 2］；

sum（a）

ans = 13 NaN 13

在向量x中刪除NaN元素，可有下列四種方法：

（1）  i=find（~isnan（x））；x=x（i）。

（2）  x=x（find（~isnan（x）））。

（3）  x=x（~isnan（x））。

（4）  x（isnan（x））=［ ］。

在矩陣X中刪除NaN所在的行，可輸入

X（any（isnan（X）‘），：）=［ ］；

經過這種預處理後的資料，可進行各種分析和統計操作。

3、迴歸和曲線擬合

對給定的資料進行擬合，可採用多項式迴歸，也可採用其它訊號形式的迴歸，其基本原理是最小二乘法，這一功能實現在MATLAB中顯得輕而易舉。

例1：設透過測量得到一組時間t與變數y的資料：

t=［0 。3 。8 1。1 1。6 2。3］；

y=［0。5 0。82 1。14 1。25 1。35 1。40］；

進行迴歸，可得到兩種不同的結果。MATLAB程式如下：

t=［0 。3 。8 1。1 1。6 2。3］’；

y=［。5 。82 1。14 1。25 1。35 1。40］‘；

X1=［ones（size（t）） t t。^2］；

a=X1\y；

X2=［ones（size（t）） exp（–t） t。*exp（–t）］；

b=X2\y；

T=［0：。1：2。5］’；

Y1=［ones（size（T）） T T。^2］*a；

Y2=［ones（size（T）） exp（-T） T。*exp（-T）］*b；

figure（1）

subplot（1，2，1）

plot（T，Y1，‘-’，t，y，‘o’），grid on

title（‘多項式迴歸’）

subplot（1，2，2）

plot（T，Y2，‘-’，t，y，‘o’），grid on

title（‘指數函式迴歸’）

例2 已知變數y與x1，x2有關，測得一組資料為

x1=［。2 。5 。6 。8 1。0 1。1 ］‘；

x2=［。1 。3 。4 。9 1。1 1。4 ］’；

y=［。17 。26 。28 。23 。27 。24］‘；

採用來擬合，則有

x1=［。2 。5 。6 。8 1。0 1。1］’；

x2=［。1 。3 。4 。9 1。1 1。4］‘；

y=［。17 。26 。28 。23 。27 。24］’；

X=［ones（size（x1）） x1 x2］；

a=X\y

a = 0。1018 0。4844 −0。2847

因此資料的擬合模型為

y=0。1018+0。4844x1−0。2487x2

4、傅立葉分析與FFT

利用MATLAB提供的FFT函式可方便地計算出訊號的傅立葉變換，從而在頻域上對訊號進行分析。

例1 ：

混合頻率訊號成分分析。有一訊號x由三種不同頻率的正弦訊號混合而成，透過得到訊號的DFT，確定出訊號的頻率及其強度關係，程式如下：

t=0：1/119：1；

x=5*sin（2*pi*20*t）+3*sin（2*pi*30*t）+sin（2*pi*45*t）；

y=fft（x）；

m=abs（y）；

f=（0：length（y） -1）‘*119/length（y）；

figure（1）

subplot（2，1，1），plot（t，x），grid on

title（’多頻率混合訊號‘）

ylabel（’Input \itx‘），xlabel（’Time ‘）

subplot（2，1，2），plot（f，m）

ylabel（’Abs。 Magnitude‘），grid on

xlabel（’Frequency （Hertz）‘）

例2 ：

訊號在傳輸過程中，由於受通道或環境影響，在接收端得到的是噪聲環境下的訊號。我們利用FFT函式對這一訊號進行傅立葉分析，從而確定訊號的頻率，程式如下：

t=0：1/199：1；

x=sin（2*pi*50*t）+1。2*randn（size（t））； ％噪聲中的訊號

y=fft（x）；

m=abs（y）；

f=（0：length（y） -1）’*199/length（y）；

figure（1）

subplot（2，1，1），plot（t，x），grid on

title（‘訊號檢測’）

ylabel（‘Input \itx’），xlabel（‘Time ’）

subplot（2，1，2），plot（f，m）

ylabel（‘Abs。 Magnitude’），grid on

xlabel（‘Frequency （Hertz）’）

例3 ：

天文學家記錄了300年來太陽黑子的活動情況，我們對這組資料進行傅立葉分析，從而得出太陽黑子的活動週期。MATLAB程式如下：

l

oad sunspot。dat

year=sunspot（：，1）；

wolfer=sunspot（：，2）；

figure（1）

subplot（2，1，1）

plot（year，wolfer）

title（‘原始資料’）

Y=fft（wolfer）；

N=length（Y）；

Y（1）=［］；

power=abs（Y（1：N/2））。^2；

nyquist=1/2；

freq=（1：N/2）/（N/2）*nyquist；

period=1。/freq；

subplot（2，1，2）

plot（period，power）

title（‘功率譜’）， grid on

axis（［0 40 0 2e7］）

各位讀者朋友，感謝您的閱讀，您若對工程應用中的數學問題感興趣，歡迎關注我，願我們一起討論和成長！！！