lab13

podsumowanie/README.md

@@ -901,11 +901,46 @@ Zagadnienia:
 
 
 ### R
+```r
+# test założeń
+shapiro.test(mtcars$wt)$p.value
+
+# pearson
+cor.test(mtcars$mpg, mtcars$wt, method = "pearson")$p.value
+cor.test(mtcars$mpg, mtcars$wt, method = "pearson")$est
+
+# kendall
+cor.test(mtcars$mpg, mtcars$wt, method = "kendall")$p.value
+cor.test(mtcars$mpg, mtcars$wt, method = "kendall")$est
+
+# spearman
+cor.test(mtcars$mpg, mtcars$wt, method = "spearman")$p.value
+cor.test(mtcars$mpg, mtcars$wt, method = "spearman")$est
+
+# inne
+chisq.test(Efekt,Metoda)
+cor(Girth,Volume)
+```
 
 
 
 ### Zagadnienia
+- Analiza zależności cech - jeśli badaniu statystycznemu podlega jednocześnie wiele cech, to jednym z podstawowych zagadnień staje się analiza zależności pomiędzy nimi. Do wykrycia zależności pomocne są odpowiednie wykresy, współczynniki mierzące jej siłę oraz testy badające jej istotność.
 
+- Tabela dwudzielcza.<br/>
+![tabela](lab13/dwudzielcza.png)
+
+- Test niezależności w tablicy dwudzielczej - po prostu wszyskie możliwe kombinacje są niezależne - na stronie jest rozpisany wzór.<br/>
+![cramer](lab13/cramer.png)
+
+- Test istotności dla współczynnika korelacji - na stronie.<br/>
+![korelacja](lab13/korelacja.png)
+
+- Korelacja Pearsona - służy do sprawdzenia czy dwie zmienne ilościowe są powiązane ze sobą związkiem liniowym. -1 = odwrotnie liniowa, 0 = brak korelacji, 1 = dodatnio liniowa.
+
+- Korelacja Tau Kendalla - statystyka będąca jedną z miar monotonicznej zależności dwóch zmiennych losowych. -1 = każda maleje przy wzroście drugiej, 0 = brak korelacji, 1 = każda ze zmiannych rośnie przy wzroście drugiej.
+
+- Korelacja rang Spearmana - ich interpretacja jest podobna do klasycznego współczynnika korelacji Pearsona, z jednym zastrzeżeniem: w odróżnieniu od współczynnika Pearsona, który mierzy liniową zależność między zmiennymi, a wszelkie inne związki traktuje jak zaburzone zależności liniowe, korelacja rangowa pokazuje dowolną monotoniczną zależność (także nieliniową).
 
 
 

testowe/.Rhistory

@@ -136,3 +136,132 @@ data.name = deparse(substitute(x)))
 class(result) <- 'htest'
 return(result)
 }
+rok <- 1995:2002
+liczba_przypadkow <- c(39.7, 38.2, 34.7, 33.1, 30.1, 28.4, 26.3, 24.7)
+data_set <- data.frame(rok = rok, liczba_przypadkow = liczba_przypadkow)
+plot(data_set, main = "Wykres rozrzutu", pch = 16)
+model <- lm(liczba_przypadkow ~ rok, data = data_set)
+model$coefficients
+plot(data_set, main = "Wykres rozrzutu", pch = 16)
+abline(model, col = "red", lwd = 2)
+coef(model)
+confint(model)
+summary(model)
+fitted(model)
+residuals(model)
+temp_rok <- data.frame(rok = seq(min(data_set$rok) - 10,
+max(data_set$rok) + 10,
+length = 100))
+pred <- stats::predict(model, temp_rok, interval = "prediction")
+plot(data_set, main = "Wykres rozrzutu", pch = 16)
+abline(model, col = "red", lwd = 2)
+lines(temp_rok$rok, pred[, 2], lty = 2, col = "red")
+lines(temp_rok$rok, pred[, 3], lty = 2, col = "red")
+temp_rok <- data.frame(rok = seq(min(data_set$rok),
+max(data_set$rok),
+length = 100))
+pred <- stats::predict(model, temp_rok, interval = "prediction")
+plot(data_set, main = "Wykres rozrzutu", pch = 16)
+abline(model, col = "red", lwd = 2)
+lines(temp_rok$rok, pred[, 2], lty = 2, col = "red")
+lines(temp_rok$rok, pred[, 3], lty = 2, col = "red")
+temp_rok <- data.frame(rok = seq(min(data_set$rok) - 10,
+max(data_set$rok) + 10,
+length = 100))
+pred <- stats::predict(model, temp_rok, interval = "prediction")
+plot(data_set, main = "Wykres rozrzutu", pch = 16)
+abline(model, col = "red", lwd = 2)
+lines(temp_rok$rok, pred[, 2], lty = 2, col = "red")
+lines(temp_rok$rok, pred[, 3], lty = 2, col = "red")
+new_rok <- data.frame(rok = 2003:2007)
+(pred_2003_2007 <- stats::predict(model, new_rok, interval = 'prediction'))
+plot(data_set, main = "Wykres rozrzutu z predykcją na lata 2003-2007", pch = 16,
+xlim = c(1995, 2007), ylim = c(10, 40))
+abline(model, col = "red", lwd = 2)
+points(2003:2007, pred_2003_2007[, 1], col = "blue", pch = 16)
+temp_rok <- data.frame(rok = seq(1994, 2008, length = 100))
+pred <- stats::predict(model, temp_rok, interval = "prediction")
+lines(temp_rok$rok, pred[, 2], lty = 2, col = "red")
+lines(temp_rok$rok, pred[, 3], lty = 2, col = "red")
+load(url("http://ls.home.amu.edu.pl/data_sets/liver_data.RData"))
+head(liver_data)
+liver_data$condition <- ifelse(liver_data$condition == "Yes", 1, 0)
+model_1 <- glm(condition ~ bilirubin + ldh, data = liver_data, family = 'binomial')
+model_1
+summary(model_1)
+step(model_1)
+exp(coef(model_1)[2])
+exp(coef(model_1)[3])
+install.packages("ROCR")
+library(ROCR)
+pred_1 <- prediction(model_1$fitted, liver_data$condition)
+plot(performance(pred_1, 'tpr', 'fpr'), main = "Model 1")
+performance(pred_1, 'auc')@y.values
+liver_data_new <- data.frame(bilirubin = c(0.9, 2.1, 3.4), ldh = c(100, 200, 300))
+(predict_glm <- stats::predict(model_1,
+liver_data_new,
+type = 'response'))
+model_1_hat <- coef(model_1)[1] +
+coef(model_1)[2] * liver_data$bilirubin +
+coef(model_1)[3] * liver_data$ldh
+model_1_temp <- seq(min(model_1_hat) - 1, max(model_1_hat) + 2.5, length.out = 100)
+condition_temp <- exp(model_1_temp) / (1 + exp(model_1_temp))
+plot(model_1_temp, condition_temp, type = "l", xlab = "X beta", ylab = "condition",
+xlim = c(-6, 9), ylim = c(-0.1, 1.1))
+points(model_1_hat, liver_data$condition, pch = 16)
+points(coef(model_1)[1] +
+coef(model_1)[2] * liver_data_new$bilirubin +
+coef(model_1)[3] * liver_data_new$ldh,
+predict_glm, pch = 16, col = "red")
+pred_1 <- prediction(model_1$fitted, liver_data$condition)
+plot(performance(pred_1, 'tpr', 'fpr'), main = "Model 1")
+performance(pred_1, 'auc')@y.values
+summary(model_1)
+model_1$coefficients
+summary(model_1)$adj.r.squared
+x1<-rexp(30,5)
+x2<-rnorm(30,2,2)
+x3<-rnorm(30,10,1)
+gen<-c(x1,x2,x3)
+# wykresy gęstości jądrowych przy różnych szerokościach okna
+par(mfrow=c(2,2))
+plot(density(gen,bw=0.1))
+plot(density(gen,bw=0.5))
+plot(density(gen,bw=3))
+plot(density(gen,bw=5))        #na trzecim
+#ZAD1
+head(USArrests)
+pairs(USArrests)
+#UrbanPop jest najsłabiej skorelowana z pozostałymi
+cor.test(USArrests$Murder,USArrests$UrbanPop, method="pearson")
+cor.test(USArrests$Rape,USArrests$UrbanPop, method="pearson")
+(pca_1 <- prcomp(~ Murder + Assault + Rape, data = USArrests, scale = TRUE))
+summary(pca_1)
+head(pca_1$x)
+cat("...")
+pca_1$rotation
+par(mfrow = c(1, 2))
+matplot(pca_1$rotation, type = 'l', lty = 1, lwd = 2,
+xlab = 'zmienne', ylab = 'ładunki', ylim = c(-0.9, 1.05),
+xaxt = "n")
+axis(1, at = 1:3, labels = rownames(pca_1$rotation))
+legend('topleft', legend = c('PC1', 'PC2', 'PC3'), ncol = 3, col = 1:3, lwd = 2)
+text(rep(1, 3), pca_1$rotation[1, ], round(pca_1$rotation[1, ], 2), pos = 4)
+text(rep(2, 3), pca_1$rotation[2, ], round(pca_1$rotation[2, ], 2), pos = 1)
+text(rep(3, 3), pca_1$rotation[3, ], round(pca_1$rotation[3, ], 2), pos = 2)
+matplot(abs(pca_1$rotation), type = 'l', lty = 1, lwd = 2,
+xlab = 'zmienne', ylab = '|ładunki|', ylim = c(0, 1.05),
+xaxt = "n")
+axis(1, at = 1:3, labels = rownames(pca_1$rotation))
+legend('topleft', legend = c('PC1', 'PC2', 'PC3'), ncol = 3, col = 1:3, lwd = 2)
+text(rep(1, 3), abs(pca_1$rotation)[1, ], abs(round(pca_1$rotation[1, ], 2)), pos = 4)
+text(rep(2, 3), abs(pca_1$rotation)[2, ], abs(round(pca_1$rotation[2, ], 2)), pos = 1)
+text(rep(3, 3), abs(pca_1$rotation)[3, ], abs(round(pca_1$rotation[3, ], 2)), pos = 2)
+plot(pca_1)
+par(mfrow = c(1, 1))
+plot(pca_1)
+pca_1$sdev^2
+mean(pca_1$sdev^2)
+biplot(pca_1)
+library(ape)
+plot(mst(dist(scale(USArrests[, -3]))), x1 = pca_1$x[, 1], x2 = pca_1$x[, 2])

zajecia14/.Rhistory

@@ -0,0 +1,143 @@
+v <- 10
+blad_i <- numeric(v)
+krok <- floor(nrow(iris) / v)
+permutacja <- sample(1:nrow(iris))
+temp <- 0
+for (i in 1:v) {
+if (i != v) {
+obs_temp <- permutacja[(temp + 1):(i * krok)]
+temp <- i * krok
+} else {
+obs_temp <- permutacja[(temp + 1):nrow(iris)]
+}
+model_lda_i <- lda(Species ~ ., data = iris[-obs_temp, ])
+pred_i <- stats::predict(model_lda_i, iris[obs_temp, ])$class
+blad_i[i] <- sum(pred_i != iris$Species[obs_temp])
+}
+sum(blad_i) / nrow(iris)
+library(MASS)
+v <- 10
+blad_i <- numeric(v)
+krok <- floor(nrow(iris) / v)
+permutacja <- sample(1:nrow(iris))
+temp <- 0
+for (i in 1:v) {
+if (i != v) {
+obs_temp <- permutacja[(temp + 1):(i * krok)]
+temp <- i * krok
+} else {
+obs_temp <- permutacja[(temp + 1):nrow(iris)]
+}
+model_lda_i <- lda(Species ~ ., data = iris[-obs_temp, ])
+pred_i <- stats::predict(model_lda_i, iris[obs_temp, ])$class
+blad_i[i] <- sum(pred_i != iris$Species[obs_temp])
+}
+sum(blad_i) / nrow(iris)
+# bootstrap
+n_boot <- 100
+temp_boot <- numeric(n_boot)
+set.seed(1234)
+for (i in 1:n_boot) {
+numery <- sample(1:nrow(iris), replace = TRUE)
+model_lda_i <- lda(Species ~ ., data = iris[numery, ])
+temp_boot[i] <- mean(stats::predict(model_lda_i, iris[-numery, ])$class != iris[-numery, ]$Species)
+}
+mean(temp_boot)
+# w pakiecie caret dla sprawdzenia
+library(caret)
+ctrl_boot <- trainControl(method = 'boot',
+number = 100,
+search = 'grid')
+ctrl_loo <- trainControl(method = 'LOOCV',
+search = 'grid')
+ctrl_10CV <- trainControl(method = "repeatedcv",
+number = 10,
+repeats = 10)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_boot)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_loo)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_10CV)
+library(caret)
+install.packages("caret")
+library(caret)
+ctrl_boot <- trainControl(method = 'boot',
+number = 100,
+search = 'grid')
+ctrl_loo <- trainControl(method = 'LOOCV',
+search = 'grid')
+ctrl_10CV <- trainControl(method = "repeatedcv",
+number = 10,
+repeats = 10)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_boot)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_loo)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_10CV)
+wina <- read.table('http://ls.home.amu.edu.pl/data_sets/wina.txt')
+head(wina)
+cat("...")
+wina$V14 <- as.factor(wina$V14)
+dim(wina)
+table(wina$V14)
+model_lda <- lda(V14 ~ V1 + V2 + V3, data = wina)
+model_lda$prior
+model_lda$means
+model_lda$scaling
+head(stats::predict(model_lda)$posterior)
+head(stats::predict(model_lda)$class)
+(conf_matrix <- table(stats::predict(model_lda)$class, wina$V14))
+(1 - sum(diag(conf_matrix)) / nrow(wina))
+(1 - sum(diag(conf_matrix)) / nrow(wina))
+pred_loo <- numeric(nrow(wina))
+for (i in 1:nrow(wina)) {
+model_lda_i <- lda(V14 ~ V1 + V2 + V3, data = wina[-i, ])
+pred_loo[i] <- stats::predict(model_lda_i, wina[i, ])$class
+}
+table(wina$V14, pred_loo)
+(1 - sum(diag(table(wina$V14, pred_loo))) / nrow(wina))
+# 10CV
+v <- 10
+blad_i <- numeric(v)
+krok <- floor(nrow(wina) / v)
+permutacja <- sample(1:nrow(wina))
+temp <- 0
+for (i in 1:v) {
+if (i != v) {
+obs_temp <- permutacja[(temp + 1):(i * krok)]
+temp <- i * krok
+} else {
+obs_temp <- permutacja[(temp + 1):nrow(wina)]
+}
+model_lda_i <- lda(V14 ~ V1 + V2 + V3, data = wina[-obs_temp, ])
+pred_i <- stats::predict(model_lda_i, wina[obs_temp, ])$class
+blad_i[i] <- sum(pred_i != wina$V14[obs_temp])
+}
+sum(blad_i) / nrow(wina)
+# bootstrap
+n_boot <- 100
+temp_boot <- numeric(n_boot)
+set.seed(1234)
+for (i in 1:n_boot) {
+numery <- sample(1:nrow(wina), replace = TRUE)
+model_lda_i <- lda(V14 ~ V1 + V2 + V3, data = wina[numery, ])
+temp_boot[i] <- mean(stats::predict(model_lda_i, wina[-numery, ])$class != wina[-numery, ]$V14)
+}
+mean(temp_boot)
+# w pakiecie caret dla sprawdzenia
+library(caret)
+ctrl_boot <- trainControl(method = 'boot',
+number = 100,
+search = 'grid')
+ctrl_loo <- trainControl(method = 'LOOCV',
+search = 'grid')
+ctrl_10CV <- trainControl(method = "repeatedcv",
+number = 10,
+repeats = 10)
+train(V14 ~ V1 + V2 + V3, data = wina, method = 'lda', trControl = ctrl_boot)
+train(V14 ~ V1 + V2 + V3, data = wina, method = 'lda', trControl = ctrl_loo)
+train(V14 ~ V1 + V2 + V3, data = wina, method = 'lda', trControl = ctrl_10CV)
+wybrane <- wina[c(20, 50, 150, 100), 1:3]
+rownames(wybrane) <- 1:4
+library(knitr)
+kable(wybrane, align = c('c', 'c', 'c'))
+new_data <- data.frame(V1 = c(13.64, 13.94, 13.08, 12.29),
+V2 = c(3.1, 1.73, 3.9, 3.17),
+V3 = c(2.56, 2.27, 2.36, 2.21))
+stats::predict(model_lda, new_data)

zajecia14/Readme.md

@@ -0,0 +1,21 @@
+# Zajęcia 14
+Klasyfikacja
+
+
+## Estymacja błędu
+Dzielimy liczbę błędów przez sumę liczby klasyfikacji.
+
+
+## Tabela krzyżowa
+Tabela przedstawiająca łączny rozkład dwóch lub większej liczby zmiennych. W przykładzie wartości w kolumnie sumują się do 100%, w wierszu już nie muszą.
+![krzyzowa](krzyzowa.png)
+
+
+## Sprawdzian krzyżowy
+Sprawdzian krzyżowy – metoda statystyczna polegająca na podziale próby statystycznej na podzbiory, a następnie przeprowadzaniu wszelkich analiz na niektórych z nich, tzw. zbiór uczący, podczas gdy pozostałe służą do potwierdzenia wiarygodności jej wyników, tzw. zbiór testowy. [link](https://pl.wikipedia.org/wiki/Sprawdzian_krzyżowy)
+
+
+## Notatki
+- Jak robię lda i dostaję wynik w którym są "coefficients of linear discriminants" to program robi tak że mnoży odpowiednią wartość współczynnika z odpowiednią wartością obserwacji i otrzymany wynik jest używany potem do algorytmu klasyfikacji.
+
+- A posteriori (łac. „z następstwa”) – w filozofii termin będący przeciwieństwem wyrażenia a priori, oznaczający tyle co: "po fakcie".

zajecia14/zadania.R

@@ -0,0 +1,146 @@
+# ZAD 1
+library(MASS)
+
+v <- 10
+blad_i <- numeric(v)
+krok <- floor(nrow(iris) / v)
+permutacja <- sample(1:nrow(iris))
+temp <- 0
+for (i in 1:v) {
+  if (i != v) {
+    obs_temp <- permutacja[(temp + 1):(i * krok)]
+    temp <- i * krok
+  } else {
+    obs_temp <- permutacja[(temp + 1):nrow(iris)]
+  }
+  model_lda_i <- lda(Species ~ ., data = iris[-obs_temp, ])
+  pred_i <- stats::predict(model_lda_i, iris[obs_temp, ])$class
+  blad_i[i] <- sum(pred_i != iris$Species[obs_temp])
+}
+sum(blad_i) / nrow(iris)
+
+
+# bootstrap
+n_boot <- 100
+temp_boot <- numeric(n_boot)
+set.seed(1234)
+for (i in 1:n_boot) {
+  numery <- sample(1:nrow(iris), replace = TRUE)
+  model_lda_i <- lda(Species ~ ., data = iris[numery, ])
+  temp_boot[i] <- mean(stats::predict(model_lda_i, iris[-numery, ])$class != iris[-numery, ]$Species)
+}
+mean(temp_boot)
+
+
+
+# w pakiecie caret dla sprawdzenia
+install.packages("caret")
+library(caret)
+ctrl_boot <- trainControl(method = 'boot',
+                          number = 100,
+                          search = 'grid')
+ctrl_loo <- trainControl(method = 'LOOCV',
+                         search = 'grid')
+ctrl_10CV <- trainControl(method = "repeatedcv",
+                          number = 10,
+                          repeats = 10)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_boot)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_loo)
+train(Species ~ ., data = iris, method = 'lda', trControl = ctrl_10CV)
+
+# ZAD 2
+wina <- read.table('http://ls.home.amu.edu.pl/data_sets/wina.txt')
+head(wina)
+cat("...")
+wina$V14 <- as.factor(wina$V14)
+
+dim(wina)
+table(wina$V14)
+
+library(MASS)
+model_lda <- lda(V14 ~ V1 + V2 + V3, data = wina)
+model_lda$prior
+model_lda$means
+model_lda$scaling
+
+head(stats::predict(model_lda)$posterior)
+head(stats::predict(model_lda)$class)
+(conf_matrix <- table(stats::predict(model_lda)$class, wina$V14))
+
+
+
+(1 - sum(diag(conf_matrix)) / nrow(wina))
+
+pred_loo <- numeric(nrow(wina))
+for (i in 1:nrow(wina)) {
+  model_lda_i <- lda(V14 ~ V1 + V2 + V3, data = wina[-i, ])
+  pred_loo[i] <- stats::predict(model_lda_i, wina[i, ])$class
+}
+table(wina$V14, pred_loo)
+(1 - sum(diag(table(wina$V14, pred_loo))) / nrow(wina))
+
+
+
+# 10CV
+v <- 10
+blad_i <- numeric(v)
+krok <- floor(nrow(wina) / v)
+permutacja <- sample(1:nrow(wina))
+temp <- 0
+for (i in 1:v) {
+  if (i != v) {
+    obs_temp <- permutacja[(temp + 1):(i * krok)]
+    temp <- i * krok
+  } else {
+    obs_temp <- permutacja[(temp + 1):nrow(wina)]
+  }
+  model_lda_i <- lda(V14 ~ V1 + V2 + V3, data = wina[-obs_temp, ])
+  pred_i <- stats::predict(model_lda_i, wina[obs_temp, ])$class
+  blad_i[i] <- sum(pred_i != wina$V14[obs_temp])
+}
+sum(blad_i) / nrow(wina)
+
+
+
+# bootstrap
+n_boot <- 100
+temp_boot <- numeric(n_boot)
+set.seed(1234)
+for (i in 1:n_boot) {
+  numery <- sample(1:nrow(wina), replace = TRUE)
+  model_lda_i <- lda(V14 ~ V1 + V2 + V3, data = wina[numery, ])
+  temp_boot[i] <- mean(stats::predict(model_lda_i, wina[-numery, ])$class != wina[-numery, ]$V14)
+}
+mean(temp_boot)
+
+
+
+# w pakiecie caret dla sprawdzenia
+library(caret)
+ctrl_boot <- trainControl(method = 'boot',
+                          number = 100,
+                          search = 'grid')
+ctrl_loo <- trainControl(method = 'LOOCV',
+                         search = 'grid')
+ctrl_10CV <- trainControl(method = "repeatedcv",
+                          number = 10,
+                          repeats = 10)
+train(V14 ~ V1 + V2 + V3, data = wina, method = 'lda', trControl = ctrl_boot)
+train(V14 ~ V1 + V2 + V3, data = wina, method = 'lda', trControl = ctrl_loo)
+train(V14 ~ V1 + V2 + V3, data = wina, method = 'lda', trControl = ctrl_10CV)
+
+
+wybrane <- wina[c(20, 50, 150, 100), 1:3]
+rownames(wybrane) <- 1:4
+library(knitr)
+kable(wybrane, align = c('c', 'c', 'c'))
+
+
+
+new_data <- data.frame(V1 = c(13.64, 13.94, 13.08, 12.29),  
+                       V2 = c(3.1, 1.73, 3.9, 3.17), 
+                       V3 = c(2.56, 2.27, 2.36, 2.21))
+stats::predict(model_lda, new_data)
+
+
+

zajecia14/zajecia14.Rproj

@@ -0,0 +1,13 @@
+Version: 1.0
+
+RestoreWorkspace: Default
+SaveWorkspace: Default
+AlwaysSaveHistory: Default
+
+EnableCodeIndexing: Yes
+UseSpacesForTab: Yes
+NumSpacesForTab: 2
+Encoding: UTF-8
+
+RnwWeave: Sweave
+LaTeX: pdfLaTeX