x <- 2*rnorm(50)
y <- x^2 - 2*x + rnorm(50)
plot(x,y)
reg1 <- lm(y ~ x)
abline(reg1$coef)
plot(x,reg1$resid)
plot(reg1$fitted,reg1$resid)
# a hibák parabolát követnek!
# regresszióba kéne négyzetes tag!
reg2 <- lm(y ~ x + I(x^2))
# I(): identitás, azért kell, mert lm()
# hasában mást jelentenek a műveleti jelek!
reg2
plot(x,y)
points(x,reg2$fitted,col="red")
lines(x,reg2$fitted,col="red")
# ha x nem növekvő sorrendben van, keszekusza...
# növekvő rendben kéne x, de ezzel együtt át kell
# rendezni a fitted vector-t is!
sort(x)
order(x)
x[order(x)]
plot(x,y)
points(x[order(x)],y[order(x)],col=2,pch=2)
lines(x[order(x)],y[order(x)],col=2)
plot(x,y)
lines(x[order(x)],reg2$fitted[order(x)],col=2)
# ha új alappontokon akarjuk ezt a becsült függvényt:
t <- seq(-6,6,0.2)
# predict függvény: illesztett regressziós függvényt
# új alappontokban képes kiszámolni
# háklis: új pontok data.frame-ben kellenek neki,
# és az oszlopok neve a becslő változók: most x
f <- predict(reg2,data.frame(x=t))
plot(x,y)
lines(t,f)

plot(x,reg2$resid)
plot(reg2$fitted,reg2$resid)
summary(reg2)
# t-teszt szerint: konstans tag nagy eséllyel 0
# új modell, konstans nélkül
reg3 <- lm(y ~ x + I(x^2) + 0)
# y = bx + cx^2 + hiba
plot(x,y)
lines(x[order(x)],reg3$fitted[order(x)])
lines(x[order(x)],reg2$fitted[order(x)],col="blue")
# grfikusan: nagyon közel vna eymáshoz a két modell/
# két függvény
# statisztikai összehasonlítás: anova
# analysis of variance
anova(reg2,reg3)
# RSS: residual sum of squares: hibák négyzetösszege
# F = F-statisztika, Pr(..)=p-value: F-teszt:
# H0: gyakorlatilag ekvivalens a két modell
# H1: lényeges a különbség
# nagy p-value: H0 -- ilyenkor az egyszerűbb modellt
# használjuk/preferáljuk; kis p-value: H1: megéri
# a komplexebb modellt használni
anova(reg1,reg2)
# anova() csak akkor van értelme, ha az egyik modell
# részmodellje a másiknak: néhány változó 
# elhagyásával kaptuk
# reg1, reg3 nem hasonlíthatók össze: mindkettőben
# van olyan változó, ami a másikban nincs.

# több változós előrejelzés:
a <- rnorm(30)
b <- rnorm(30)
c <- 2*a - b + rnorm(30)
reg4 <- lm(c ~ a + b)
reg4
# itt is lehet szorzat, magasabb rendű tag(ok), ...
reg5 <- lm(c ~ a + b + I(a*b) + I(a^2))