290 likes | 404 Views
Statick é modely-aproximácia polynómom. Všeobecný tvar polynómu:. Statické modely. p = polyfit(x,y,n) % x,y merané udaje [p,S] = polyfit(x,y,n) % S chyba odhadu [p,S,mu] = polyfit(x,y,n) [p,S,mu] = polyfit(x,y,n) %odhad koeficientov polynomu. % Príklad1
E N D
Statické modely-aproximácia polynómom Všeobecný tvar polynómu:
Statické modely • p = polyfit(x,y,n) % x,y merané udaje • [p,S] = polyfit(x,y,n)% S chyba odhadu • [p,S,mu] = polyfit(x,y,n) • [p,S,mu] = polyfit(x,y,n) %odhad koeficientov polynomu
% Príklad1 x = (0: 0.1: 2.5)';% Vstupné údaje: y = erf(x); % erf(x) = 2/sqrt(pi) * integral from 0 to x of exp(-t^2) dt. plot(x,y) grid pause p = polyfit(x,y,6) % najdenie koef. polynomu 6.radu pause f = polyval(p,x); % vypocet modelovanych hodnot table = [x y f y-f] % vytvorenie tabulky hodnot mer.+odchylky pause plot(x,y,'o',x,f,'+')
Priklad 2 Staicky model 3.radu x = [1 2 3 4 5]'; %merana nez. premenna y = [5.5 43.1 128 290.7 498.4]'; %merany vystup p = polyfit(x,y,3) ym = polyval(p,x); % vypocet vyst. z modelu table = [x y ym y-ym] % vytvorenie tabulky pause plot(x,y,'o',x,ym,'+') pause x2 = [1:.1:5]';% generovane nove x y2 = polyval(p,x2); %vypocet vystupu s novymi x plot(x,y,'o',x2,y2, '-') grid
Regresná analýza Merané údaje: t = [0 .3 .8 1.1 1.6 2.3]'; y = [0.5 0.82 1.14 1.25 1.35 1.40]'; plot(t,y,'o'), grid on Tvar regresnej krivky (yM)
X = [ones(size(t)) t t.^2] X = 1.0000 0 0 1.0000 0.3000 0.0900 1.0000 0.8000 0.6400 1.0000 1.1000 1.2100 1.0000 1.6000 2.5600 1.0000 2.3000 5.2900 The solution is found with the backslash operator. a = X\y
a = X\y a = 0.5318 0.9191 - 0.2387 The second-order polynomial model of the data is therefore Now evaluate the model at regularly spaced points and overlay the original data in a plot. T = (0:0.1:2.5)'; Y = [ones(size(T)) T T.^2]*a; plot(T,Y,'-',t,y,'o'), grid on
Linear-in-the-Parameters Regression Tvar regresnej krivky : X = [ones(size(t)) exp(-t) t.*exp(-t)]; a = X\y a = 1.3974 - 0.8988 0.4097
X = [ones(size(t)) exp(-t) t.*exp(-t)] X = 1.0000 1.0000 0 1.0000 0.7408 0.2222 1.0000 0.4493 0.3595 1.0000 0.3329 0.3662 1.0000 0.2019 0.3230 1.0000 0.1003 0.2306
Mnohonásobná regresia Výstupná veličina y=f(x1,x2,...,xn)je funkciou viacerých vstupných premenných. Príklad: Pozorovania (merania) x1 = [.2 .5 .6 .8 1.0 1.1]';% prvý vstup x2 = [.1 .3 .4 .9 1.1 1.4]';% druhý vstup y = [.17 .26 .28 .23 .27 .24]';% výstup Model y=a0+a1x1+a2x2 Neznáme koeficientya0, a1,a2(metóda najmenších štvorcov) , formovaním regresnej matice X, X = [ones(size(x1)) x1 x2]; a = X\y
a = 0.1018 0.4844 -0.2847 Hľadaný model: y=0.1018+0.4844x1-0.2847x2 Y = X*a; MaxErr = max(abs(Y - y)) MaxErr = 0.0038
Priklad 4 (nahrada parabolou) t = [0 .3 .8 1.1 1.6 2.3]'; % vstup y = [0.5 0.82 1.14 1.25 1.35 1.40]'; % vystup plot(t,y,'o'), grid on pause X = [ones(size(t)) t t.^2] % regresna matica a = X\y % optimalne parametre, preurc. system rovnic T = (0:0.1:2.5)'; % test. modelu na novych udajoch Y = [ones(size(T)) T T.^2]*a; plot(T,Y,'-',t,y,'o') grid on
Nahrada exponencialou X = [ones(size(t)) exp(- t) t.*exp(- t)]; a = X\y pause T = (0:0.1:2.5)'; Y = [ones(size(T)) exp(- T) T.*exp(- T)]*a; plot(T,Y, '+',t,y,'o') grid on
Mnohonásobná regresia x1 = [.2 .5 .6 .8 1.0 1.1]'; % prvy vstup x2 = [.1 .3 .4 .9 1.1 1.4]'; % druhy vstup y = [.17 .26 .28 .23 .27 .24]'; % vystup pause X = [ones(size(x1)) x1 x2]; a = X\y % optimalne riesenie-parametre Y = X*a; Maxchyba = max(abs(Y - y))
a = 1; b = [1/4 1/4 1/4 1/4]; % 4 spriemerovane hodnoty load count.dat % nacitanie udajov x = count(:,1); % nefiltrovana premenna y = filter(b,a,x); % filtrovaná premenná t = 1:length(x); plot(t,x,'-.',t,y,'-') grid on legend(' Povodne udaje', ' Filtrovane udaje',2)
x = [4 3 7 -9 1 0 0 0]' ; y = fft(x) load sunspot.dat % analyza aktivity sl.skvrn year = sunspot(:,1); wolfer = sunspot(:,2); plot(year,wolfer) title('Sunspot Data') 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; plot(freq,power), grid on xlabel('cycles/year') title('Periodogram')
load sunspot.dat % analyza aktivity sl.skvrn year = sunspot(:,1); wolfer = sunspot(:,2); plot(year,wolfer) title('Sunspot Data'),display('Pozri graf') pause 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; plot(freq,power), grid on xlabel('cycles/year') title('Periodogram')
period = 1./freq; plot(period,power), axis([0 40 0 2e7]), grid on ylabel('Power') xlabel('Period(Years/Cycle)') pause [mp index] = max(power); period(index)
t = 0:1/100:10-1/100; x = sin(2*pi*15*t) + sin(2*pi*40*t); y = fft(x); m = abs(y); p = unwrap(angle(y)); f = (0:length(y)-1)'*100/length(y); subplot(2,1,1), plot(f,m), ylabel('Abs. Magnitude'), grid on subplot(2,1,2), plot(f,p*180/pi) ylabel('Phase [Degrees]'), grid on xlabel('Frequency [Hertz]')