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.:

Logistička regresija

Izlaz funkcije uvijek je između 0 i 1. Provjerite sliku ispod

Logistička regresija

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

Logistička regresija

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
  1. 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:

Provjerite kontinuirane varijable

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
  1. 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]]

Provjerite varijable faktora

## ## [[2]]

Provjerite varijable faktora

## ## [[3]]

Provjerite varijable faktora

## ## [[4]]

Provjerite varijable faktora

## ## [[5]]

Provjerite varijable faktora

## ## [[6]]

Provjerite varijable faktora

Napomena: Koristite sljedeći gumb za navigaciju do sljedećeg grafikona

Provjerite varijable faktora

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:

Sumarna statistika

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:

Sumarna statistika

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:

Sumarna statistika

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:

Sumarna statistika

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:

Nelinearnost

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:

Korelacija

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.

Matrica zabune

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

Matrica zabune

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

  1. 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 Preciznost protiv prisjećanja rezultat na temelju preciznosti i prisjećanja. The Preciznost protiv prisjećanja je harmonijska sredina ove dvije metrike, što znači da daje veću težinu nižim vrijednostima.

Preciznost protiv prisjećanja

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:

ROC krivulja

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”)