GLM u R: Generalizirani linearni model s primjerom
Što je logistička regresija?
Logistička regresija se koristi za predviđanje klase, tj. vjerojatnosti. Logistička regresija može točno predvidjeti binarni ishod.
Zamislite da želite predvidjeti hoće li zajam biti odbijen/prihvaćen na temelju mnogih atributa. Logistička regresija je oblika 0/1. y = 0 ako je posudba odbijena, y = 1 ako je prihvaćena.
Logistički regresijski model razlikuje se od linearnog regresijskog modela na dva načina.
- Prije svega, logistička regresija prihvaća samo dihotomni (binarni) input kao zavisnu varijablu (tj. vektor od 0 i 1).
- Drugo, ishod se mjeri sljedećom funkcijom vjerojatnosne veze tzv sigmoidni zbog svog oblika slova S.:
Izlaz funkcije uvijek je između 0 i 1. Provjerite sliku ispod
Sigmoidna funkcija vraća vrijednosti od 0 do 1. Za zadatak klasifikacije potreban nam je diskretni izlaz od 0 ili 1.
Za pretvaranje kontinuiranog protoka u diskretnu vrijednost, možemo postaviti ograničenje odluke na 0.5. Sve vrijednosti iznad ovog praga klasificiraju se kao 1
Kako stvoriti generalizirani model linije (GLM)
Koristimo se odrasla osoba skup podataka za ilustraciju logističke regresije. "Odrasli" je izvrstan skup podataka za zadatak klasifikacije. Cilj je predvidjeti hoće li godišnji prihod u dolarima pojedinca premašiti 50.000. Skup podataka sadrži 46,033 XNUMX opažanja i deset značajki:
- starost: starost pojedinca. Numerički
- obrazovanje: Obrazovna razina pojedinca. Faktor.
- bračni.status: Maritalni status pojedinca. Čimbenik tj. nikad u braku, u braku-građanski supružnik, …
- spol: Spol pojedinca. Faktor tj. Muško ili Žensko
- prihod: Target varijabla. Prihod iznad ili ispod 50K. Faktor tj. >50K, <=50K
između ostalog
library(dplyr) data_adult <-read.csv("https://raw.githubusercontent.com/guru99-edu/R-Programming/master/adult.csv") glimpse(data_adult)
Izlaz:
Observations: 48,842 Variables: 10 $ x <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,... $ age <int> 25, 38, 28, 44, 18, 34, 29, 63, 24, 55, 65, 36, 26... $ workclass <fctr> Private, Private, Local-gov, Private, ?, Private,... $ education <fctr> 11th, HS-grad, Assoc-acdm, Some-college, Some-col... $ educational.num <int> 7, 9, 12, 10, 10, 6, 9, 15, 10, 4, 9, 13, 9, 9, 9,... $ marital.status <fctr> Never-married, Married-civ-spouse, Married-civ-sp... $ race <fctr> Black, White, White, Black, White, White, Black, ... $ gender <fctr> Male, Male, Male, Male, Female, Male, Male, Male,... $ hours.per.week <int> 40, 50, 40, 40, 30, 30, 40, 32, 40, 10, 40, 40, 39... $ income <fctr> <=50K, <=50K, >50K, >50K, <=50K, <=50K, <=50K, >5...
Postupit ćemo na sljedeći način:
- Korak 1: Provjerite kontinuirane varijable
- Korak 2: Provjerite faktorske varijable
- Korak 3: Inženjering značajki
- Korak 4: Zbirna statistika
- Korak 5: Skup za obuku/testiranje
- Korak 6: Izgradite model
- Korak 7: Procijenite izvedbu modela
- korak 8: Poboljšajte model
Vaš zadatak je predvidjeti koja će osoba imati prihod veći od 50K.
U ovom vodiču detaljno će biti opisan svaki korak za izvođenje analize na stvarnom skupu podataka.
Korak 1) Provjerite kontinuirane varijable
U prvom koraku možete vidjeti distribuciju kontinuiranih varijabli.
continuous <-select_if(data_adult, is.numeric) summary(continuous)
Objašnjenje koda
- kontinuirano <- select_if(data_adult, is.numeric): Koristite funkciju select_if() iz biblioteke dplyr za odabir samo numeričkih stupaca
- sažetak (kontinuirano): Ispis sažete statistike
Izlaz:
## X age educational.num hours.per.week ## Min. : 1 Min. :17.00 Min. : 1.00 Min. : 1.00 ## 1st Qu.:11509 1st Qu.:28.00 1st Qu.: 9.00 1st Qu.:40.00 ## Median :23017 Median :37.00 Median :10.00 Median :40.00 ## Mean :23017 Mean :38.56 Mean :10.13 Mean :40.95 ## 3rd Qu.:34525 3rd Qu.:47.00 3rd Qu.:13.00 3rd Qu.:45.00 ## Max. :46033 Max. :90.00 Max. :16.00 Max. :99.00
Iz gornje tablice možete vidjeti da podaci imaju potpuno različite skale i sati.po.tjednima imaju velike odstupanja (tj. pogledajte zadnji kvartil i maksimalnu vrijednost).
S tim se možete nositi u dva koraka:
- 1: Nacrtajte distribuciju sati po tjednu
- 2: Standardizirajte kontinuirane varijable
- Nacrtajte distribuciju
Pogledajmo pobliže distribuciju sati po tjednu
# Histogram with kernel density curve library(ggplot2) ggplot(continuous, aes(x = hours.per.week)) + geom_density(alpha = .2, fill = "#FF6666")
Izlaz:
Varijabla ima puno outliera i nije dobro definiranu distribuciju. Ovaj problem možete djelomično riješiti brisanjem gornjih 0.01 posto sati tjedno.
Osnovna sintaksa kvantila:
quantile(variable, percentile) arguments: -variable: Select the variable in the data frame to compute the percentile -percentile: Can be a single value between 0 and 1 or multiple value. If multiple, use this format: `c(A,B,C, ...) - `A`,`B`,`C` and `...` are all integer from 0 to 1.
Izračunavamo gornja 2 posto percentila
top_one_percent <- quantile(data_adult$hours.per.week, .99) top_one_percent
Objašnjenje koda
- kvantil(data_adult$hours.per.week, .99): Izračunajte vrijednost 99 posto radnog vremena
Izlaz:
## 99% ## 80
98 posto stanovništva radi manje od 80 sati tjedno.
Možete ispustiti opažanja iznad ovog praga. Koristite filtar iz dplyr knjižnica.
data_adult_drop <-data_adult %>% filter(hours.per.week<top_one_percent) dim(data_adult_drop)
Izlaz:
## [1] 45537 10
- Standardizirajte kontinuirane varijable
Možete standardizirati svaki stupac kako biste poboljšali izvedbu jer vaši podaci nemaju istu ljestvicu. Možete koristiti funkciju mutate_if iz biblioteke dplyr. Osnovna sintaksa je:
mutate_if(df, condition, funs(function)) arguments: -`df`: Data frame used to compute the function - `condition`: Statement used. Do not use parenthesis - funs(function): Return the function to apply. Do not use parenthesis for the function
Brojčane stupce možete standardizirati na sljedeći način:
data_adult_rescale <- data_adult_drop % > % mutate_if(is.numeric, funs(as.numeric(scale(.)))) head(data_adult_rescale)
Objašnjenje koda
- mutate_if(is.numeric, funs(scale)): uvjet je samo numerički stupac, a funkcija je scale
Izlaz:
## X age workclass education educational.num ## 1 -1.732680 -1.02325949 Private 11th -1.22106443 ## 2 -1.732605 -0.03969284 Private HS-grad -0.43998868 ## 3 -1.732530 -0.79628257 Local-gov Assoc-acdm 0.73162494 ## 4 -1.732455 0.41426100 Private Some-college -0.04945081 ## 5 -1.732379 -0.34232873 Private 10th -1.61160231 ## 6 -1.732304 1.85178149 Self-emp-not-inc Prof-school 1.90323857 ## marital.status race gender hours.per.week income ## 1 Never-married Black Male -0.03995944 <=50K ## 2 Married-civ-spouse White Male 0.86863037 <=50K ## 3 Married-civ-spouse White Male -0.03995944 >50K ## 4 Married-civ-spouse Black Male -0.03995944 >50K ## 5 Never-married White Male -0.94854924 <=50K ## 6 Married-civ-spouse White Male -0.76683128 >50K
Korak 2) Provjerite faktorske varijable
Ovaj korak ima dva cilja:
- Provjerite razinu u svakom stupcu kategorije
- Definirajte nove razine
Ovaj korak ćemo podijeliti u tri dijela:
- Odaberite kategoričke stupce
- Pohranite trakasti grafikon svakog stupca na popis
- Ispišite grafikone
Stupce faktora možemo odabrati pomoću donjeg koda:
# Select categorical column factor <- data.frame(select_if(data_adult_rescale, is.factor)) ncol(factor)
Objašnjenje koda
- data.frame(select_if(data_adult, is.factor)): Stupce faktora pohranjujemo u faktor u tipu okvira podataka. Biblioteka ggplot2 zahtijeva objekt okvira podataka.
Izlaz:
## [1] 6
Skup podataka sadrži 6 kategoričkih varijabli
Drugi korak je vještiji. Želite iscrtati trakasti grafikon za svaki stupac u faktoru okvira podataka. Pogodnije je automatizirati proces, posebno u situaciji kada postoji mnogo stupaca.
library(ggplot2) # Create graph for each column graph <- lapply(names(factor), function(x) ggplot(factor, aes(get(x))) + geom_bar() + theme(axis.text.x = element_text(angle = 90)))
Objašnjenje koda
- lapply(): Koristite funkciju lapply() za prosljeđivanje funkcije u svim stupcima skupa podataka. Izlaz pohranjujete na popis
- funkcija(x): funkcija će biti obrađena za svaki x. Ovdje je x stupac
- ggplot(faktor, aes(get(x))) + geom_bar()+ theme(axis.text.x = element_text(angle = 90)): Napravite trakasti dijagram za svaki x element. Napomena, da vratite x kao stupac, morate ga uključiti unutar get()
Posljednji korak je relativno lak. Želite ispisati 6 grafikona.
# Print the graph graph
Izlaz:
## [[1]]
## ## [[2]]
## ## [[3]]
## ## [[4]]
## ## [[5]]
## ## [[6]]
Napomena: Koristite sljedeći gumb za navigaciju do sljedećeg grafikona
Korak 3) Inženjering značajki
Preinačiti obrazovanje
Iz gornjeg grafikona možete vidjeti da varijabla obrazovanje ima 16 razina. To je značajno, a neke razine imaju relativno mali broj opažanja. Ako želite poboljšati količinu informacija koju možete dobiti iz ove varijable, možete je preinačiti na višu razinu. Naime, stvarate veće grupe sa sličnim stupnjem obrazovanja. Na primjer, niska razina obrazovanja pretvorit će se u napuštanje škole. Više razine obrazovanja bit će promijenjene u master.
Evo detalja:
Stara razina | Nova razina |
---|---|
Predškolski | maloljetnik koji se ispisao iz škole |
10th | maloljetnik koji se ispisao iz škole |
11th | maloljetnik koji se ispisao iz škole |
12th | maloljetnik koji se ispisao iz škole |
1.-4 | maloljetnik koji se ispisao iz škole |
5th-6th | maloljetnik koji se ispisao iz škole |
7th-8th | maloljetnik koji se ispisao iz škole |
9th | maloljetnik koji se ispisao iz škole |
HS-Grad | HighGrad |
Neki fakulteti | Zajednica |
izv.-acdm | Zajednica |
izv. voc | Zajednica |
Celibat | Celibat |
Masters | Masters |
Prof-škola | Masters |
Doktorat | Dr. |
recast_data <- data_adult_rescale % > % select(-X) % > % mutate(education = factor(ifelse(education == "Preschool" | education == "10th" | education == "11th" | education == "12th" | education == "1st-4th" | education == "5th-6th" | education == "7th-8th" | education == "9th", "dropout", ifelse(education == "HS-grad", "HighGrad", ifelse(education == "Some-college" | education == "Assoc-acdm" | education == "Assoc-voc", "Community", ifelse(education == "Bachelors", "Bachelors", ifelse(education == "Masters" | education == "Prof-school", "Master", "PhD")))))))
Objašnjenje koda
- Koristimo glagol mutirati iz biblioteke dplyr. Izjavom ifelse mijenjamo vrijednosti obrazovanja
U donjoj tablici kreirate sažetak statistike da vidite, u prosjeku, koliko je godina obrazovanja (z-vrijednost) potrebno da se postigne diploma, magisterij ili doktorat.
recast_data % > % group_by(education) % > % summarize(average_educ_year = mean(educational.num), count = n()) % > % arrange(average_educ_year)
Izlaz:
## # A tibble: 6 x 3 ## education average_educ_year count ## <fctr> <dbl> <int> ## 1 dropout -1.76147258 5712 ## 2 HighGrad -0.43998868 14803 ## 3 Community 0.09561361 13407 ## 4 Bachelors 1.12216282 7720 ## 5 Master 1.60337381 3338 ## 6 PhD 2.29377644 557
Preoblikovati Marital-status
Također je moguće stvoriti niže razine za bračni status. U sljedećem kodu mijenjate razinu na sljedeći način:
Stara razina | Nova razina |
---|---|
Nikad oženjen | Neoženjen |
Oženjen-supružnik-odsutan | Neoženjen |
Oženjen-AF-bračni drug | Oženjen |
Oženjen-građanka-supružnik | |
Odvojen | Odvojen |
Rastavljen | |
udovice | Udovica |
# Change level marry recast_data <- recast_data % > % mutate(marital.status = factor(ifelse(marital.status == "Never-married" | marital.status == "Married-spouse-absent", "Not_married", ifelse(marital.status == "Married-AF-spouse" | marital.status == "Married-civ-spouse", "Married", ifelse(marital.status == "Separated" | marital.status == "Divorced", "Separated", "Widow")))))
Možete provjeriti broj pojedinaca unutar svake grupe.
table(recast_data$marital.status)
Izlaz:
## ## Married Not_married Separated Widow ## 21165 15359 7727 1286
Korak 4) Sažeta statistika
Vrijeme je da provjerimo neke statistike o našim ciljnim varijablama. Na donjem grafikonu brojite postotak pojedinaca koji zarađuju više od 50 tisuća s obzirom na njihov spol.
# Plot gender income ggplot(recast_data, aes(x = gender, fill = income)) + geom_bar(position = "fill") + theme_classic()
Izlaz:
Zatim provjerite utječe li porijeklo pojedinca na njihovu zaradu.
# Plot origin income ggplot(recast_data, aes(x = race, fill = income)) + geom_bar(position = "fill") + theme_classic() + theme(axis.text.x = element_text(angle = 90))
Izlaz:
Broj sati rada prema spolu.
# box plot gender working time ggplot(recast_data, aes(x = gender, y = hours.per.week)) + geom_boxplot() + stat_summary(fun.y = mean, geom = "point", size = 3, color = "steelblue") + theme_classic()
Izlaz:
Okvirni dijagram potvrđuje da raspodjela radnog vremena odgovara različitim skupinama. U okvirnom dijagramu, oba spola nemaju homogena opažanja.
Gustoću tjednog radnog vremena možete provjeriti po vrsti obrazovanja. Distribucije imaju mnogo različitih odabira. To se vjerojatno može objasniti vrstom ugovora u SAD-u.
# Plot distribution working time by education ggplot(recast_data, aes(x = hours.per.week)) + geom_density(aes(color = education), alpha = 0.5) + theme_classic()
Objašnjenje koda
- ggplot(recast_data, aes( x= hours.per.week)): Dijagram gustoće zahtijeva samo jednu varijablu
- geom_density(aes(color = education), alpha =0.5): Geometrijski objekt za kontrolu gustoće
Izlaz:
Kako biste potvrdili svoje misli, možete izvesti jednosmjerni ANOVA test:
anova <- aov(hours.per.week~education, recast_data) summary(anova)
Izlaz:
## Df Sum Sq Mean Sq F value Pr(>F) ## education 5 1552 310.31 321.2 <2e-16 *** ## Residuals 45531 43984 0.97 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ANOVA test potvrđuje razliku u prosjeku između skupina.
Nelinearnost
Prije nego što pokrenete model, možete vidjeti je li broj radnih sati povezan s dobi.
library(ggplot2) ggplot(recast_data, aes(x = age, y = hours.per.week)) + geom_point(aes(color = income), size = 0.5) + stat_smooth(method = 'lm', formula = y~poly(x, 2), se = TRUE, aes(color = income)) + theme_classic()
Objašnjenje koda
- ggplot(recast_data, aes(x = age, y = sati.po.tjednu)): Postavite estetiku grafikona
- geom_point(aes(color= prihod), veličina =0.5): Konstruirajte točkasti dijagram
- stat_smooth(): Dodajte liniju trenda sa sljedećim argumentima:
- method='lm': iscrtajte prilagođenu vrijednost ako je Linearna regresija
- formula = y~poly(x,2): Uklopi polinomsku regresiju
- se = TRUE: Dodajte standardnu pogrešku
- aes(boja= prihod): Razdvojite model prema prihodu
Izlaz:
Ukratko, možete testirati uvjete interakcije u modelu kako biste uočili učinak nelinearnosti između tjednog radnog vremena i drugih značajki. Važno je otkriti pod kojim se uvjetima radno vrijeme razlikuje.
Korelacija
Sljedeća provjera je vizualizacija korelacije između varijabli. Vrstu razine faktora pretvarate u numeričku tako da možete iscrtati toplinsku kartu koja sadrži koeficijent korelacije izračunat Spearmanovom metodom.
library(GGally) # Convert data to numeric corr <- data.frame(lapply(recast_data, as.integer)) # Plot the graphggcorr(corr, method = c("pairwise", "spearman"), nbreaks = 6, hjust = 0.8, label = TRUE, label_size = 3, color = "grey50")
Objašnjenje koda
- data.frame(lapply(recast_data,as.integer)): Pretvori podatke u numeričke
- ggcorr() crta toplinsku kartu sa sljedećim argumentima:
- metoda: Metoda za izračunavanje korelacije
- nbreaks = 6: Broj prekida
- hjust = 0.8: Kontrolni položaj naziva varijable u dijagramu
- label = TRUE: Dodajte oznake u središte prozora
- label_size = 3: Oznake veličine
- boja = “siva50”): Boja naljepnice
Izlaz:
Korak 5) Skup za obuku/testiranje
Svaki pod nadzorom stroj za učenje Zadatak zahtijeva dijeljenje podataka između skupa vlakova i testnog skupa. Možete upotrijebiti "funkciju" koju ste izradili u drugim vodičima za nadzirano učenje za izradu skupa za treniranje/testiranje.
set.seed(1234) create_train_test <- function(data, size = 0.8, train = TRUE) { n_row = nrow(data) total_row = size * n_row train_sample <- 1: total_row if (train == TRUE) { return (data[train_sample, ]) } else { return (data[-train_sample, ]) } } data_train <- create_train_test(recast_data, 0.8, train = TRUE) data_test <- create_train_test(recast_data, 0.8, train = FALSE) dim(data_train)
Izlaz:
## [1] 36429 9
dim(data_test)
Izlaz:
## [1] 9108 9
Korak 6) Izgradite model
Da biste vidjeli kako algoritam radi, koristite paket glm(). The Generalizirani linearni model je kolekcija modela. Osnovna sintaksa je:
glm(formula, data=data, family=linkfunction() Argument: - formula: Equation used to fit the model- data: dataset used - Family: - binomial: (link = "logit") - gaussian: (link = "identity") - Gamma: (link = "inverse") - inverse.gaussian: (link = "1/mu^2") - poisson: (link = "log") - quasi: (link = "identity", variance = "constant") - quasibinomial: (link = "logit") - quasipoisson: (link = "log")
Spremni ste procijeniti logistički model kako biste podijelili razinu prihoda između skupa značajki.
formula <- income~. logit <- glm(formula, data = data_train, family = 'binomial') summary(logit)
Objašnjenje koda
- formula <- prihod ~ .: Stvorite model koji odgovara
- logit <- glm(formula, data = data_train, family = 'binomial'): Uklopi logistički model (family = 'binomial') s podacima data_train.
- summary(logit): Ispis sažetka modela
Izlaz:
## ## Call: ## glm(formula = formula, family = "binomial", data = data_train) ## ## Deviance Residuals: ## Min 1Q Median 3Q Max ## -2.6456 -0.5858 -0.2609 -0.0651 3.1982 ## ## Coefficients: ## Estimate Std. Error z value Pr(>|z|) ## (Intercept) 0.07882 0.21726 0.363 0.71675 ## age 0.41119 0.01857 22.146 < 2e-16 *** ## workclassLocal-gov -0.64018 0.09396 -6.813 9.54e-12 *** ## workclassPrivate -0.53542 0.07886 -6.789 1.13e-11 *** ## workclassSelf-emp-inc -0.07733 0.10350 -0.747 0.45499 ## workclassSelf-emp-not-inc -1.09052 0.09140 -11.931 < 2e-16 *** ## workclassState-gov -0.80562 0.10617 -7.588 3.25e-14 *** ## workclassWithout-pay -1.09765 0.86787 -1.265 0.20596 ## educationCommunity -0.44436 0.08267 -5.375 7.66e-08 *** ## educationHighGrad -0.67613 0.11827 -5.717 1.08e-08 *** ## educationMaster 0.35651 0.06780 5.258 1.46e-07 *** ## educationPhD 0.46995 0.15772 2.980 0.00289 ** ## educationdropout -1.04974 0.21280 -4.933 8.10e-07 *** ## educational.num 0.56908 0.07063 8.057 7.84e-16 *** ## marital.statusNot_married -2.50346 0.05113 -48.966 < 2e-16 *** ## marital.statusSeparated -2.16177 0.05425 -39.846 < 2e-16 *** ## marital.statusWidow -2.22707 0.12522 -17.785 < 2e-16 *** ## raceAsian-Pac-Islander 0.08359 0.20344 0.411 0.68117 ## raceBlack 0.07188 0.19330 0.372 0.71001 ## raceOther 0.01370 0.27695 0.049 0.96054 ## raceWhite 0.34830 0.18441 1.889 0.05894 . ## genderMale 0.08596 0.04289 2.004 0.04506 * ## hours.per.week 0.41942 0.01748 23.998 < 2e-16 *** ## ---## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## (Dispersion parameter for binomial family taken to be 1) ## ## Null deviance: 40601 on 36428 degrees of freedom ## Residual deviance: 27041 on 36406 degrees of freedom ## AIC: 27087 ## ## Number of Fisher Scoring iterations: 6
Sažetak našeg modela otkriva zanimljive podatke. Učinkovitost logističke regresije procjenjuje se pomoću specifičnih ključnih metrika.
- AIC (Akaike informacijski kriteriji): Ovo je ekvivalent R2 u logističkoj regresiji. Mjeri prikladnost kada se kazna primijeni na broj parametara. Manji AIC vrijednosti pokazuju da je model bliži istini.
- Nulta devijacija: odgovara modelu samo s presjekom. Stupanj slobode je n-1. Možemo je protumačiti kao vrijednost hi-kvadrat (prilagođena vrijednost koja se razlikuje od testiranja hipoteze stvarne vrijednosti).
- Rezidualna devijacija: Model sa svim varijablama. Također se tumači kao testiranje hi-kvadrat hipoteze.
- Broj ponavljanja Fisherova bodovanja: Broj ponavljanja prije konvergiranja.
Izlaz funkcije glm() pohranjuje se na popis. Kod u nastavku prikazuje sve stavke dostupne u logit varijabli koju smo konstruirali za procjenu logističke regresije.
# Popis je jako dug, ispišite samo prva tri elementa
lapply(logit, class)[1:3]
Izlaz:
## $coefficients ## [1] "numeric" ## ## $residuals ## [1] "numeric" ## ## $fitted.values ## [1] "numeric"
Svaka se vrijednost može izdvojiti znakom $ nakon kojeg slijedi naziv metrike. Na primjer, pohranili ste model kao logit. Za izdvajanje AIC kriterija koristite:
logit$aic
Izlaz:
## [1] 27086.65
Korak 7) Procijenite izvedbu modela
Matrica zabune
The matrica zabune je bolji izbor za procjenu izvedbe klasifikacije u usporedbi s različitim mjernim podacima koje ste vidjeli prije. Opća ideja je brojati koliko su puta True instance klasificirane kao False.
Da biste izračunali matricu zabune, prvo morate imati skup predviđanja kako bi se mogla usporediti sa stvarnim ciljevima.
predict <- predict(logit, data_test, type = 'response') # confusion matrix table_mat <- table(data_test$income, predict > 0.5) table_mat
Objašnjenje koda
- predict(logit,data_test, type = 'response'): Izračunajte predviđanje na skupu testova. Postavite type = 'response' za izračunavanje vjerojatnosti odgovora.
- table(data_test$income, predict > 0.5): Izračunajte matricu zabune. predict > 0.5 znači da vraća 1 ako su predviđene vjerojatnosti iznad 0.5, inače 0.
Izlaz:
## ## FALSE TRUE ## <=50K 6310 495 ## >50K 1074 1229
Svaki redak u matrici zabune predstavlja stvarni cilj, dok svaki stupac predstavlja predviđeni cilj. Prvi redak ove matrice uzima u obzir dohodak niži od 50 tisuća (klasa lažnih): 6241 je ispravno klasificirano kao pojedinci s dohotkom nižim od 50 tisuća (Istinski negativan), dok je preostali pogrešno klasificiran kao iznad 50 tisuća (Lažno pozitivno). U drugom redu su primanja iznad 50k, pozitivna klasa je 1229 (Prava pozitiva), dok Istinski negativan bio 1074.
Možete izračunati model točnost zbrajanjem stvarnog pozitivnog + istinskog negativnog ukupnog promatranja
accuracy_Test <- sum(diag(table_mat)) / sum(table_mat) accuracy_Test
Objašnjenje koda
- sum(diag(table_mat)): Zbroj dijagonale
- sum(table_mat): Zbroj matrice.
Izlaz:
## [1] 0.8277339
Čini se da model pati od jednog problema, precjenjuje broj lažno negativnih rezultata. Ovo se zove paradoks testa točnosti. Naveli smo da je točnost omjer točnih predviđanja prema ukupnom broju slučajeva. Možemo imati relativno visoku točnost, ali beskoristan model. To se događa kada postoji dominantna klasa. Ako pogledate unatrag na matricu zabune, možete vidjeti da je većina slučajeva klasificirana kao stvarno negativna. Zamislite sada, model je klasificirao sve klase kao negativne (tj. niže od 50k). Imali biste točnost od 75 posto (6718/6718+2257). Vaš model ima bolju izvedbu, ali muči se razlikovati pravo pozitivno od pravog negativnog.
U takvoj situaciji poželjno je imati sažetiju metriku. Možemo pogledati:
- Preciznost=TP/(TP+FP)
- Opoziv=TP/(TP+FN)
Preciznost protiv prisjećanja
Preciznost gleda na točnost pozitivnog predviđanja. Podsjetiti je omjer pozitivnih instanci koje je klasifikator ispravno otkrio;
Možete konstruirati dvije funkcije za izračunavanje ove dvije metrike
- Konstruirajte preciznost
precision <- function(matrix) { # True positive tp <- matrix[2, 2] # false positive fp <- matrix[1, 2] return (tp / (tp + fp)) }
Objašnjenje koda
- mat[1,1]: vraća prvu ćeliju prvog stupca podatkovnog okvira, tj. pravi pozitivan
- prostirka[1,2]; Vrati prvu ćeliju drugog stupca podatkovnog okvira, tj. lažno pozitivno
recall <- function(matrix) { # true positive tp <- matrix[2, 2]# false positive fn <- matrix[2, 1] return (tp / (tp + fn)) }
Objašnjenje koda
- mat[1,1]: vraća prvu ćeliju prvog stupca podatkovnog okvira, tj. pravi pozitivan
- prostirka[2,1]; Vrati drugu ćeliju prvog stupca podatkovnog okvira, tj. lažno negativno
Možete testirati svoje funkcije
prec <- precision(table_mat) prec rec <- recall(table_mat) rec
Izlaz:
## [1] 0.712877 ## [2] 0.5336518
Kada model kaže da se radi o osobi iznad 50 tisuća, točan je u samo 54 posto slučajeva, a može potraživati osobe iznad 50 tisuća u 72 posto slučajeva.
Možete stvoriti rezultat na temelju preciznosti i prisjećanja. The
je harmonijska sredina ove dvije metrike, što znači da daje veću težinu nižim vrijednostima.
f1 <- 2 * ((prec * rec) / (prec + rec)) f1
Izlaz:
## [1] 0.6103799
Kompromis preciznosti u odnosu na opoziv
Nemoguće je imati i visoku preciznost i visoko pamćenje.
Ako povećamo preciznost, točna će se osoba bolje predvidjeti, ali bismo propustili puno njih (slabije pamćenje). U nekim situacijama više volimo veću preciznost od opoziva. Postoji konkavan odnos između preciznosti i prisjećanja.
- Zamislite, morate predvidjeti ima li pacijent bolest. Želite biti što precizniji.
- Ako trebate otkriti potencijalne prevarante na ulici putem prepoznavanja lica, bilo bi bolje uhvatiti mnogo ljudi označenih kao prevaranti iako je preciznost niska. Policija će moći pustiti osobu koja nije prijevaru.
ROC krivulja
The Prijamnik Operating Karakterističan krivulja je još jedan uobičajeni alat koji se koristi s binarnom klasifikacijom. Vrlo je slična krivulji preciznost/prisjećanje, ali umjesto iscrtavanja preciznosti u odnosu na prisjećanje, ROC krivulja pokazuje pravu pozitivnu stopu (tj. prisjećanje) naspram lažno pozitivne stope. Lažno pozitivna stopa je omjer negativnih slučajeva koji su netočno klasificirani kao pozitivni. Jednaka je jedan minus prava negativna stopa. Prava negativna stopa također se naziva specifičnost. Stoga je prikazana ROC krivulja osjetljivost (prisjećanje) u odnosu na 1-specifičnost
Za iscrtavanje ROC krivulje moramo instalirati biblioteku koja se zove RORC. Možemo pronaći u kondi knjižnica. Možete upisati kod:
conda install -cr r-rocr –da
Možemo iscrtati ROC s funkcijama prediction() i performance().
library(ROCR) ROCRpred <- prediction(predict, data_test$income) ROCRperf <- performance(ROCRpred, 'tpr', 'fpr') plot(ROCRperf, colorize = TRUE, text.adj = c(-0.2, 1.7))
Objašnjenje koda
- predviđanje(predviđanje, data_test$income): ROCR biblioteka treba stvoriti objekt predviđanja za transformaciju ulaznih podataka
- performance(ROCRpred, 'tpr','fpr'): Vrati dvije kombinacije za izradu na grafikonu. Ovdje su konstruirani tpr i fpr. Za precizno iscrtavanje i opoziv zajedno koristite "prec", "rec".
Izlaz:
Korak 8) Poboljšajte model
Možete pokušati dodati nelinearnost modelu s interakcijom između
- dob i sati.po.tjednu
- spol i sati.po.tjednu.
Morate upotrijebiti test rezultata da biste usporedili oba modela
formula_2 <- income~age: hours.per.week + gender: hours.per.week + . logit_2 <- glm(formula_2, data = data_train, family = 'binomial') predict_2 <- predict(logit_2, data_test, type = 'response') table_mat_2 <- table(data_test$income, predict_2 > 0.5) precision_2 <- precision(table_mat_2) recall_2 <- recall(table_mat_2) f1_2 <- 2 * ((precision_2 * recall_2) / (precision_2 + recall_2)) f1_2
Izlaz:
## [1] 0.6109181
Ocjena je nešto viša od prethodne. Možete nastaviti raditi na podacima pokušavajući nadmašiti rezultat.
rezime
Možemo sažeti funkciju za treniranje logističke regresije u tablici u nastavku:
Paket | Cilj | funkcija | Argument |
---|---|---|---|
- | Stvorite skup podataka o vlaku/testiranju | create_train_set() | podaci, veličina, vlak |
glm | Uvježbajte generalizirani linearni model | glm() | formula, podaci, obitelj* |
glm | Sažmite model | Sažetak() | ugrađeni model |
baza | Napravite predviđanje | predvidjeti() | prilagođeni model, skup podataka, tip = 'odgovor' |
baza | Napravite matricu zabune | stol() | y, predvidi() |
baza | Stvorite ocjenu točnosti | zbroj(dijag(tablica())/zbroj(tablica() | |
ROCR | Stvorite ROC: Korak 1 Stvorite predviđanje | predviđanje() | predvidjeti(), y |
ROCR | Stvorite ROC: 2. korak Stvorite izvedbu | izvođenje() | predviđanje(), 'tpr', 'fpr' |
ROCR | Stvorite ROC: korak 3 Iscrtajte grafikon | zemljište() | izvođenje() |
Drugi GLM vrste modela su:
– binom: (link = “logit”)
– gaussov: (link = “identitet”)
– Gama: (link = “obrnuto”)
– inverse.gaussian: (link = “1/mu^2”)
– poisson: (link = “log”)
– kvazi: (link = “identitet”, varijanca = “konstanta”)
– kvazibinom: (link = “logit”)
– quasipoisson: (link = “log”)