# Lab: Classification Methods ## The Stock Market Data ### library(ISLR2) ?Smarket head(Smarket) names(Smarket) dim(Smarket) summary(Smarket) ### pairs(Smarket) # Direction == "up" akko Today>0 cor(Smarket[, -9]) # velika korelacija izmedu Volume i Year; mala korelacija Today i Lag_k boxplot(Volume~Year, data=Smarket) attach(Smarket) ### train <- (Year < 2005) Smarket.2005 <- Smarket[!train, ] dim(Smarket.2005) # 252 povrata u 2005 n=1250-252 ## velicina skupa za trening Direction.2005 <- Direction[!train] ############################ logisticka regresija ################################### # Neka je Y=Direction te pretpostavljamo model # P(Y="Up" | Lag1,...., Lag5, Volume) = h(b0 + b1*Lag1 + ... + b5*Lag5 + b6*Volume), za # h(z) = exp(z)/(1+exp(z)) glm.fits <- glm( Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + Volume, data = Smarket, family = binomial, subset = train ) summary(glm.fits) ## procijenjeni parametri (b0_hat,..., b6_hat) su dakle coef(glm.fits) # p-vrijednost u zadnjem stupcu odnosi se na testiranje H_0: bj=0 pri cemu je testna # statistika slicna kao kod lin. regresije (tzv. Waldova statistika), te rezultate slicno intepretiramo # -> niti jedna se zasebno ne cini znacajna (najmanja p-vrijednost za Lag1) # glm.fits0 <- glm( # Direction ~ 1, # data = Smarket, family = binomial, subset = train # ) # # anova(glm.fits0, glm.fits, test="chisq") # ### velika p-vrijednost -- ne mozemo odbaciti nultu hipotezu da je model samo sa slobodnim clanom # ### (dakle, bez kovarijata) dovoljan. ## procjena testne greske glm.probs <- predict(glm.fits, Smarket.2005, type = "response") glm.pred <- rep("Down", 252) glm.pred[glm.probs > .5] <- "Up" table(glm.pred, Direction.2005) # kontingencijska tablica za vektor (glm.pred, Direction.2005) # Direction.2005 # glm.pred Down Up # Down 77 97 # Up 34 44 ### procjena testne greske P(Y_hat != Y) (dakle, uz 0-1 funkciju gubitka) je mean(glm.pred != Direction.2005) ## 0.5198413 = (34+97)/252 --> losije od slucajnog pogadanja! ## probajmo: glm.fits <- glm(Direction ~ Lag1 + Lag2, data = Smarket, family = binomial, subset = train) glm.probs <- predict(glm.fits, Smarket.2005, type = "response") glm.pred <- rep("Down", 252) glm.pred[glm.probs > .5] <- "Up" table(glm.pred, Direction.2005) # Direction.2005 # glm.pred Down Up # Down 35 35 # Up 76 106 mean(glm.pred != Direction.2005) # [1] 0.4404762 --> izgleda puno bolje, ali # kada bi stavili f_hat(x):= "up" za sve x, dobili bi istz testnu gresku (to je vjerojatno slucajnost) mean(Direction.2005 != "Up") # [1] 0.4404762 # Ipak, uocimo da je u nasem modelu procjena za uvjetnu gresku P(Y != Y_hat | Y_hat="Up") 76/(76+106) # [1] 0.4175824 ## Primjer predikcije h=function(z) exp(z)/(1+exp(z)) # nas model je sada # P(Y="Up" | Lag1, Lag2) = h(b0 + b1*Lag1 + b2*Lag2) b=coef(glm.fits) # procjena za P(Y="Up" | Lag1=1.2, Lag2=1.1) je h(b[1] + b[2]*1.2 + b[3]*1.1) # 0.4791462 # ili preko predict: predict(glm.fits, newdata = data.frame(Lag1 = 1.2, Lag2 = 1.1), type = "response" ) # 0.4791462 ######################################### LDA ############################### library(MASS) lda.fit <- lda(Direction ~ Lag1 + Lag2, data = Smarket, subset = train) # sintaksa je ista kao kod glm() funkcije, samo bez opcije "family" lda.fit # dakle, ako je 0="Down", 1="Up", # 1. pi0_hat=0.491984, pi1_hat=0.508016 # 2. mi0_hat=(0.04279022, 0.03389409) approx E(X | Y=0), mi1_hat=(-0.03954635 -0.03132544) # 3. Coefficients of linear discriminants -- o tome cemo malo kasnije lda.pred <- predict(lda.fit, Smarket.2005) names(lda.pred) # lda.pred$posterior daje P(Y= "Down" | X=x) i P(Y= "Up" | X=x) lda.class <- lda.pred$class table(lda.class, Direction.2005) mean(lda.class != Direction.2005) # [1] 0.4404762 -- isto kao kod log. regresije jer daju potpuno iste predikcije na testnom skupu: table(lda.class, glm.pred) # ### Coefficients of linear discriminants: # w=as.vector(lda.fit$scaling) ## <-- vektor koeficijenata (ili smjer) prve disktriminacijske koordinate (jedina jer K=1) # pi=as.vector(lda.fit$prior) # mi_0=as.vector(lda.fit$means[1,]) # mi_1=as.vector(lda.fit$means[2,]) # # lda.pred <- predict(lda.fit, Smarket.2005)$class # # daje procijenjene klase na testnom skupu # # ### pokazimo da je # ### f_hat(x) = argmax_{k=0,1} delta_k(x) = argmin_{k=0,1} [ ||w*x - w*mi_k_hat ||/2 - log(pi_k_hat) ] # delta_0=function(x){abs(w%*%x - w%*%mi_0)^2/2 - log(pi[1])} # delta_1=function(x){abs(w%*%x - w%*%mi_1)^2/2 - log(pi[2])} # f_hat=function(x){ # if (delta_0(x) Ovo nije dobra procjena za testnu gresku! ### Opcenito, skup za odabir parametra kompleksnosti i procjenu greske ne smije biti isti! ### Koristimo dodatni skup za tzv. validaciju.