Regresja liniowa i dopasowanie innych funkcji

Clear[x,y,α,β,γ]

Konstrukcja próbki

Na różne sposoby konstruuję losową próbke "danych" o rozmiarze n.

n=40;

Próbka 1

Rozkład jednostajny w polu {-1,1}x{-1,1}

X=RandomReal[{-1,1},{n,2}];
p1=ListPlot[X]

Próbka 2

Rozkład jednostajny w otoczeniu o szerokosci ϵ prostej y=x

a=2; b=3;
ϵ=0.1;

XX:=Module[{},xx=RandomReal[{a,b},1];yy=xx+RandomReal[{-ϵ,ϵ},1]; {xx,yy}//Flatten]
X=Table[XX,{i,1,n}];

(Uwaga: Nie objaśniam polecenia Module w tym miejscu!)

p1=ListPlot[X]

Próbka 3

Rozkład jednostajny w otoczeniu o szerokosci ϵ sinusoidy  y= Sin[x],  w przedziale x ∈ {a, b}

a=2; b=3;
ϵ=0.1;

XX:=Module[{},xx=RandomReal[{a,b},1];yy=Sin[xx]+RandomReal[{-ϵ,ϵ},1]; {xx,yy}//Flatten]
X=Table[XX,{i,1,n}];

(Uwaga: Nie objaśniam polecenia Module w tym miejscu!)

p1=ListPlot[X]

Regresja liniowa

Funkcja LinearModelFit

X

y=Normal[LinearModelFit[X,x,x]]

p2=Plot[y,{x,a,b}];

Show[{p1,p2}]

Odstępstwa od prawa liniowego

LinearModelFit[X,x,x]["FitResiduals"]
ListPlot[%, Filling->Axis]

Metoda najmniejszych kwadratów (suma kwadratów)

X

TableForm[X,TableDirections->Column,TableHeadings->{{"pomiar"},{"x","y"}}]

Szukam funkcji  y=α x +β  dla której suma kwadratów odchyleń punktów pomiarowych od prostej jest minimalna.

𝕩=Transpose[X][[1]];
𝕪=Transpose[X][[2]];

Lista błędów odchylenia:

ϵ=𝕪-α 𝕩-β //Simplify

F[α_,β_]=(∑_(i=1)^n)((ϵ[[i]])^2)                                              //Simplify
F[α,β]==(∑_(i=1)^n)(((𝕪[[i]]-α 𝕩[[i]]-β))^2)                    //Simplify

F[α,β]==ϵ.ϵ                                                            //Simplify

Metoda najmniejszych kwadratów (minimalizacja)

Minimalizuję błąd

warunek={(∂_α)F[α,β]==0,(∂_β)F[α,β]==0};

wspolczynniki=Solve[warunek,{α,β}]//Flatten

Sprawdzam zgodność moich obliczeń z wynikiem procedury LinearModelFit

y-(α x+β /.wspolczynniki)//Chop

Metoda najmniejszych kwadratów (niepewność pomiarowa)

Niepewnosci standardowe {α ,β} 
wg. B. Kamys: Statystyczne Metody Opracowania Pomiarów (I)
http://users.uj.edu.pl/~ufkamys/BK/smop1N_p.pdf

Show[{p1,p2}]

W=n  𝕩.𝕩-((Total[𝕩])^2)

σ[ϵ]=StandardDeviation[ϵ/.wspolczynniki]

σ[a]=σ[ϵ](Sqrt[(n/W)])

σ[b]=σ[ϵ](Sqrt[((𝕩.𝕩)/W)])

Dopasowanie (dowolnych) funkcji

Funkcja FindFit

Konstruuję wielomian n-tego stopnia

s=2;
A=Table["α"<>ToString[k],{k,0,s}]//ToExpression
P[x_]=(∑_(k=0)^s)A[[k+1]] (x^k)

FindFit[X,P[x], A,x];
ff=P[x]/.%

p3=Plot[ ff,{x,a,b}];

Show[{p1,p3}]