Vhodná polynomiální regrese v R

. (Můžete nahlásit problém o obsahu na této stránce zde)chcete sdílet svůj obsah na R-bloggerů? klikněte zde, pokud máte blog, nebo zde, pokud nemáte.

lineární vztah mezi dvěma proměnnými x a y je jedním z nejběžnějších, nejúčinnějších a nejjednodušších předpokladů, které je třeba učinit, když se snažíte zjistit jejich vztah. Někdy však, skutečný základní vztah je složitější než to, a to je, když polynomiální regrese přichází na pomoc.

podívejme se na příklad z ekonomie: Předpokládejme, že byste chtěli koupit určité množství q určitého produktu. Pokud je jednotková cena p, pak byste zaplatili celkovou částku y. Toto je typický příklad lineárního vztahu. Celková cena a množství jsou přímo úměrné. Abychom to vykreslili, napsali bychom něco takového:

p <- 0.5q <- seq(0,100,1)y <- p*qplot(q,y,type='l',col='red',main='Linear relationship')

děj bude vypadat takto:
linear-relationship

Toto je dobrá aproximace skutečného vztahu mezi y a q, ale při nákupu a prodeji bychom mohli zvážit některé další relevantní informace, například: Nákup značného množství je pravděpodobné, že můžeme požádat a získat slevu, nebo kupovat více a více určitého zboží, které bychom mohli tlačit cenu nahoru.
to může vést ke scénáři, jako je tento, kdy celkové náklady již nejsou lineární funkcí množství:

y <- 450 + p*(q-10)^3plot(q,y,type='l',col='navy',main='Nonlinear relationship',lwd=3)

Rplot02
pomocí polynomiální regrese můžeme do dat vložit modely řádu n > 1 a pokusit se modelovat nelineární vztahy.

jak přizpůsobit polynomiální regresi

nejprve nezapomeňte použít set.seed(n) při generování pseudonáhodných čísel. Tímto způsobem generátor náhodných čísel generuje vždy stejná čísla.

set.seed(20)

prediktor (q). Použijte seq pro generování rovnoměrně rozložených sekvencí rychle

q <- seq(from=0, to=20, by=0.1)

hodnota pro predikci (y):

y <- 500 + 0.4 * (q-10)^3

nějaký šum je generován a přidán do reálného signálu (y):

noise <- rnorm(length(q), mean=10, sd=80)noisy.y <- y + noise

děj hlučného signálu:

plot(q,noisy.y,col='deepskyblue4',xlab='q',main='Observed data')lines(q,y,col='firebrick1',lwd=3)

toto je graf našich simulovaných pozorovaných dat. Simulované datové body jsou modré tečky, zatímco červená čára je signál (signál je technický termín, který se často používá k označení obecného trendu, který nás zajímá).
Rplot03
náš model by měl být něco takového: y = a * q + b * q2 + c * q3 + náklady

pojďme to přizpůsobit pomocí R. Při montáži polynomů můžete buď použít

model <- lm(noisy.y ~ poly(q,3))

nebo

model <- lm(noisy.y ~ x + I(X^2) + I(X^3))

Mějte však na paměti, že q, I(q^2) a I(q^3) budou korelovány a korelované proměnné mohou způsobit problémy. Použití poly() vám umožní vyhnout se tím, že vytvoří ortogonální polynomy, proto budu používat první možnost.

summary(model)Call:lm(formula = noisy.y ~ poly(q, 3))Residuals: Min 1Q Median 3Q Max -212.326 -51.186 4.276 61.485 165.960 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 513.615 5.602 91.69 <2e-16 ***poly(q, 3)1 2075.899 79.422 26.14 <2e-16 ***poly(q, 3)2 -108.004 79.422 -1.36 0.175 poly(q, 3)3 864.025 79.422 10.88 <2e-16 ***---Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Residual standard error: 79.42 on 197 degrees of freedomMultiple R-squared: 0.8031,Adjusted R-squared: 0.8001 F-statistic: 267.8 on 3 and 197 DF, p-value: 0 

pomocí funkce confint() můžeme získat intervaly spolehlivosti parametrů našeho modelu.
intervaly spolehlivosti parametrů modelu:

confint(model, level=0.95) 2.5 % 97.5 %(Intercept) 502.5676 524.66261poly(q, 3)1 1919.2739 2232.52494poly(q, 3)2 -264.6292 48.62188poly(q, 3)3 707.3999 1020.65097

spiknutí namontovaných vs zbytků. Pokud je model vhodný, neměl by se na zbytkovém grafu zobrazit žádný jasný vzor

plot(fitted(model),residuals(model))

Rplot04
celkově se model jeví jako dobrý fit, jak ukazuje r na druhou 0,8. Koeficienty výrazu prvního a třetího řádu jsou statisticky významné, jak jsme očekávali. Nyní můžeme pomocí funkce predict() získat nastavené hodnoty a intervaly spolehlivosti, abychom vše vykreslili proti našim datům.

předpokládané hodnoty a intervaly spolehlivosti:

predicted.intervals <- predict(model,data.frame(x=q),interval='confidence', level=0.99)

přidat řádky do stávajícího grafu:

lines(q,predicted.intervals,col='green',lwd=3)lines(q,predicted.intervals,col='black',lwd=1)lines(q,predicted.intervals,col='black',lwd=1)

přidat legendu:

legend("bottomright",c("Observ.","Signal","Predicted"), col=c("deepskyblue4","red","green"), lwd=3)

zde je spiknutí:
Rplot05
vidíme, že náš model odvedl slušnou práci při montáži dat, a proto s ním můžeme být spokojeni.

varování: polynomy jsou mocné nástroje, ale mohou selhat: v tomto případě jsme věděli, že původní signál byl generován pomocí polynomu třetího stupně, ale při analýze reálných dat o něm obvykle víme jen málo, a proto musíme být opatrní, protože použití polynomů vysokého řádu (n > 4) může vést k nadměrné montáži. K nadměrné montáži dochází, když váš model místo signálu zvedá šum: i když se váš model zlepšuje a lépe přizpůsobuje stávající data, to může být špatné, když se snažíte předvídat nová data a vést k zavádějícím výsledkům.

podstatu s úplným kódem pro tento příklad naleznete zde.

Děkujeme, že jste si přečetli tento příspěvek, zanechte komentář níže, pokud máte nějaké dotazy.

Napsat komentář

Vaše e-mailová adresa nebude zveřejněna.