GLM i R: Generaliseret lineær model med eksempel
Evelyn Clarke
Hvad er logistisk regression?
Logistisk regression bruges til at forudsige en klasse, dvs. en sandsynlighed. Logistisk regression kan forudsige et binært resultat nøjagtigt.
Forestil dig, at du vil forudsige, om et lån nægtes/accepteres baseret på mange egenskaber. Den logistiske regression er af formen 0/1. y = 0 hvis et lån afvises, y = 1 hvis det accepteres.
En logistisk regressionsmodel adskiller sig fra lineær regressionsmodel på to måder.
- Først og fremmest accepterer den logistiske regression kun dikotomisk (binær) input som en afhængig variabel (dvs. en vektor på 0 og 1).
- For det andet måles resultatet af følgende probabilistiske linkfunktion kaldet sigmoid på grund af sin S-form.:
Funktionens output er altid mellem 0 og 1. Tjek billede nedenfor
Sigmoid-funktionen returnerer værdier fra 0 til 1. Til klassificeringsopgaven har vi brug for et diskret output på 0 eller 1.
For at konvertere et kontinuerligt flow til diskret værdi kan vi sætte en beslutningsgrænse til 0.5. Alle værdier over denne tærskel er klassificeret som 1
Sådan opretter du Generalized Liner Model (GLM)
Lad os bruge voksen datasæt til at illustrere logistisk regression. "Voksen" er et fantastisk datasæt til klassificeringsopgaven. Målet er at forudsige, om en persons årlige indkomst i dollar vil overstige 50.000. Datasættet indeholder 46,033 observationer og ti funktioner:
- alder: individets alder. Numerisk
- uddannelse: Den enkeltes uddannelsesniveau. Faktor.
- ægteskabelig.status: Mariden enkeltes status. Faktor dvs. aldrig gift, gift-civ-ægtefælle, …
- køn: Individets køn. Faktor, dvs. mand eller kvinde
- indkomst: Target variabel. Indkomst over eller under 50K. Faktor dvs. >50K, <=50K
blandt andre
library(dplyr)
data_adult <-read.csv("https://raw.githubusercontent.com/guru99-edu/R-Programming/master/adult.csv")
glimpse(data_adult)
Output:
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...
Vi fortsætter som følger:
- Trin 1: Tjek kontinuerte variabler
- Trin 2: Tjek faktorvariabler
- Trin 3: Feature engineering
- Trin 4: Sammenfattende statistik
- Trin 5: Træn/testsæt
- Trin 6: Byg modellen
- Trin 7: Vurder modellens ydeevne
- trin 8: Forbedre modellen
Din opgave er at forudsige, hvilken person der vil have en omsætning på mere end 50K.
I denne vejledning vil hvert trin blive detaljeret for at udføre en analyse på et rigtigt datasæt.
Trin 1) Tjek kontinuerte variable
I det første trin kan du se fordelingen af de kontinuerte variable.
continuous <-select_if(data_adult, is.numeric) summary(continuous)
Kode Forklaring
- kontinuerlig <- select_if(data_adult, is.numeric): Brug funktionen select_if() fra dplyr-biblioteket til kun at vælge de numeriske kolonner
- summary (kontinuerlig): Udskriv oversigtsstatistikken
Output:
## 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
Fra ovenstående tabel kan du se, at dataene har helt forskellige skalaer og timer.per.uger har store outliers (.dvs. se på sidste kvartil og maksimumværdi).
Du kan håndtere det ved at følge to trin:
- 1: Tegn fordelingen af timer.per.uge
- 2: Standardiser de kontinuerte variable
- Tegn fordelingen
Lad os se nærmere på fordelingen af timer.per.uge
# Histogram with kernel density curve
library(ggplot2)
ggplot(continuous, aes(x = hours.per.week)) +
geom_density(alpha = .2, fill = "#FF6666")
Output:
Variablen har masser af outliers og ikke veldefineret fordeling. Du kan delvist løse dette problem ved at slette de øverste 0.01 procent af timerne om ugen.
Grundlæggende syntaks for kvantil:
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.
Vi beregner den øverste 2 procent percentil
top_one_percent <- quantile(data_adult$hours.per.week, .99) top_one_percent
Kode Forklaring
- quantile(data_adult$hours.per.week, .99): Beregn værdien af de 99 procent af arbejdstiden
Output:
## 99% ## 80
98 procent af befolkningen arbejder under 80 timer om ugen.
Du kan slippe observationerne over denne tærskel. Du bruger filteret fra dplyr bibliotek.
data_adult_drop <-data_adult %>% filter(hours.per.week<top_one_percent) dim(data_adult_drop)
Output:
## [1] 45537 10
- Standardiser de kontinuerte variable
Du kan standardisere hver kolonne for at forbedre ydeevnen, fordi dine data ikke har samme skala. Du kan bruge funktionen mutate_if fra dplyr-biblioteket. Den grundlæggende syntaks er:
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
Du kan standardisere de numeriske kolonner som følger:
data_adult_rescale <- data_adult_drop % > % mutate_if(is.numeric, funs(as.numeric(scale(.)))) head(data_adult_rescale)
Kode Forklaring
- mutate_if(is.numeric, funs(scale)): Betingelsen er kun numerisk kolonne, og funktionen er skala
Output:
## 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
Trin 2) Tjek faktorvariabler
Dette trin har to mål:
- Tjek niveauet i hver kategorisk kolonne
- Definer nye niveauer
Vi vil opdele dette trin i tre dele:
- Vælg de kategoriske kolonner
- Gem søjlediagrammet for hver kolonne på en liste
- Udskriv graferne
Vi kan vælge faktorkolonnerne med nedenstående kode:
# Select categorical column factor <- data.frame(select_if(data_adult_rescale, is.factor)) ncol(factor)
Kode Forklaring
- data.frame(select_if(data_adult, is.factor)): Vi gemmer faktorkolonnerne i faktor i en datarammetype. Biblioteket ggplot2 kræver et datarammeobjekt.
Output:
## [1] 6
Datasættet indeholder 6 kategoriske variable
Det andet trin er mere dygtig. Du vil plotte et søjlediagram for hver kolonne i datarammefaktoren. Det er mere bekvemt at automatisere processen, især i situationer, hvor der er mange kolonner.
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)))
Kode Forklaring
- lapply(): Brug funktionen lapply() til at sende en funktion i alle kolonnerne i datasættet. Du gemmer outputtet på en liste
- funktion(x): Funktionen vil blive behandlet for hver x. Her er x kolonnerne
- ggplot(faktor, aes(get(x))) + geom_bar()+ theme(axis.text.x = element_text(angle = 90)): Opret et søjlediagram for hvert x-element. Bemærk, for at returnere x som en kolonne, skal du inkludere det i get()
Det sidste trin er relativt nemt. Du vil udskrive de 6 grafer.
# Print the graph graph
Output:
## [[1]]
## ## [[2]]
## ## [[3]]
## ## [[4]]
## ## [[5]]
## ## [[6]]
Bemærk: Brug næste knap til at navigere til næste graf
Trin 3) Feature engineering
Omarbejdet uddannelse
Af grafen ovenfor kan du se, at den variable uddannelse har 16 niveauer. Dette er betydeligt, og nogle niveauer har et relativt lavt antal observationer. Hvis du ønsker at forbedre mængden af information, du kan få fra denne variabel, kan du omforme den til et højere niveau. Man opretter nemlig større grupper med tilsvarende uddannelsesniveau. Eksempelvis vil lavt uddannelsesniveau blive konverteret til frafald. Højere uddannelsesniveauer ændres til master.
Her er detaljen:
| Gammelt niveau | Nyt niveau |
|---|---|
| Børnehave | droppe ud |
| 10. | Droppe ud |
| 11. | Droppe ud |
| 12. | Droppe ud |
| 1.-4 | Droppe ud |
| 5th-6th | Droppe ud |
| 7th-8th | Droppe ud |
| 9. | Droppe ud |
| HS-Grad | HighGrad |
| Et eller andet universitet | Community |
| Assoc-acdm | Community |
| Assoc-voc | Community |
| bachelorer | bachelorer |
| Masters | Masters |
| Prof-skole | Masters |
| Doktorgrad | PhD |
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")))))))
Kode Forklaring
- Vi bruger verbet mutere fra dplyr bibliotek. Vi ændrer uddannelsens værdier med udsagnet ifelse
I nedenstående tabel laver du en opsummerende statistik for i gennemsnit at se, hvor mange års uddannelse (z-værdi), der skal til for at nå Bachelor, Master eller PhD.
recast_data % > % group_by(education) % > % summarize(average_educ_year = mean(educational.num), count = n()) % > % arrange(average_educ_year)
Output:
## # 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
Omarbejdning Marital-status
Det er også muligt at lave lavere niveauer for civilstanden. I følgende kode ændrer du niveauet som følger:
| Gammelt niveau | Nyt niveau |
|---|---|
| Aldrig gift | Ikke gift |
| Gift-ægtefælle-fraværende | Ikke gift |
| Gift-AF-ægtefælle | Gift |
| Gift-civ-ægtefælle | |
| Separeret | Separeret |
| Fraskilt | |
| Enker | Enke |
# 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")))))
Du kan kontrollere antallet af personer inden for hver gruppe.
table(recast_data$marital.status)
Output:
## ## Married Not_married Separated Widow ## 21165 15359 7727 1286
Trin 4) Sammenfattende statistik
Det er på tide at tjekke nogle statistikker om vores målvariabler. I grafen nedenfor tæller du procentdelen af personer, der tjener mere end 50, givet deres køn.
# Plot gender income
ggplot(recast_data, aes(x = gender, fill = income)) +
geom_bar(position = "fill") +
theme_classic()
Output:
Dernæst skal du kontrollere, om individets oprindelse påvirker deres indtjening.
# 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))
Output:
Antallet af arbejdstimer fordelt på køn.
# 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()
Output:
Boksplottet bekræfter, at fordelingen af arbejdstid passer til forskellige grupper. I boksplotten har begge køn ikke homogene observationer.
Du kan tjekke tætheden af den ugentlige arbejdstid efter uddannelsestype. Fordelingerne har mange forskellige valg. Det kan formentlig forklares med kontrakttypen i USA.
# Plot distribution working time by education
ggplot(recast_data, aes(x = hours.per.week)) +
geom_density(aes(color = education), alpha = 0.5) +
theme_classic()
Kode Forklaring
- ggplot(recast_data, aes( x= hours.per.week)): Et tæthedsplot kræver kun én variabel
- geom_density(aes(color = uddannelse), alpha =0.5): Det geometriske objekt til at styre tætheden
Output:
For at bekræfte dine tanker kan du udføre en envejs ANOVA test:
anova <- aov(hours.per.week~education, recast_data) summary(anova)
Output:
## 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-testen bekræfter forskellen i gennemsnit mellem grupperne.
Ikke-linearitet
Inden du kører modellen, kan du se, om antallet af arbejdstimer er relateret til alder.
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()
Kode Forklaring
- ggplot(recast_data, aes(x = alder, y = timer.per.uge)): Indstil grafens æstetik
- geom_point(aes(farve=indkomst), størrelse =0.5): Konstruer prikplottet
- stat_smooth(): Tilføj trendlinjen med følgende argumenter:
- method='lm': Plot den tilpassede værdi, hvis lineær regression
- formel = y~poly(x,2): Tilpas en polynomiel regression
- se = TRUE: Tilføj standardfejlen
- aes(farve=indkomst): Bryd modellen efter indkomst
Output:
I en nøddeskal kan du teste interaktionsbegreber i modellen for at opfange den ikke-linearitetseffekt mellem den ugentlige arbejdstid og andre funktioner. Det er vigtigt at opdage, under hvilke forhold arbejdstiden adskiller sig.
Korrelation
Næste kontrol er at visualisere sammenhængen mellem variablerne. Du konverterer faktorniveautypen til numerisk, så du kan plotte et varmekort, der indeholder korrelationskoefficienten beregnet med Spearman-metoden.
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")
Kode Forklaring
- data.frame(lapply(recast_data,as.integer)): Konverter data til numeriske
- ggcorr() plot varmekortet med følgende argumenter:
- metode: Metode til at beregne korrelationen
- nbreaks = 6: Antal pauser
- hjust = 0.8: Kontrolposition for variabelnavnet i plottet
- label = TRUE: Tilføj labels i midten af vinduerne
- label_size = 3: Størrelsesetiketter
- farve = "grå50"): Farve på etiketten
Output:
Trin 5) Træn/testsæt
Enhver overvåget machine learning opgave kræver at opdele data mellem et togsæt og et testsæt. Du kan bruge den "funktion", du oprettede i de andre superviserede læringsvejledninger, til at oprette et tog-/testsæt.
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)
Output:
## [1] 36429 9
dim(data_test)
Output:
## [1] 9108 9
Trin 6) Byg modellen
For at se hvordan algoritmen fungerer, bruger du glm()-pakken. Det Generaliseret lineær model er en samling af modeller. Den grundlæggende syntaks er:
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")
Du er klar til at estimere den logistiske model for at opdele indkomstniveauet mellem et sæt funktioner.
formula <- income~. logit <- glm(formula, data = data_train, family = 'binomial') summary(logit)
Kode Forklaring
- formel <- indkomst ~ .: Opret modellen, så den passer
- logit <- glm(formel, data = data_train, family = 'binomial'): Tilpas en logistisk model (family = 'binomial') med data_train-dataene.
- summary(logit): Udskriv resuméet af modellen
Output:
## ## 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
Sammenfatningen af vores model afslører interessante oplysninger. Ydeevnen af en logistisk regression evalueres med specifikke nøglemålinger.
- AIC (Akaike Information Criteria): Dette svarer til R2 i logistisk regression. Den måler pasformen, når der pålægges en straf på antallet af parametre. Mindre AIC værdier indikerer, at modellen er tættere på sandheden.
- Nul afvigelse: Passer kun til modellen med intercept. Frihedsgraden er n-1. Vi kan fortolke det som en chi-kvadratværdi (fitted værdi forskellig fra den faktiske værdi hypotesetest).
- Residual Deviance: Model med alle variablerne. Det tolkes også som en chi-kvadrat hypotesetestning.
- Antal Fisher Scoring iterationer: Antal iterationer før konvergering.
Outputtet af glm()-funktionen er gemt i en liste. Koden nedenfor viser alle de tilgængelige elementer i den logit-variabel, vi konstruerede for at evaluere den logistiske regression.
# Listen er meget lang, udskriv kun de tre første elementer
lapply(logit, class)[1:3]
Output:
## $coefficients ## [1] "numeric" ## ## $residuals ## [1] "numeric" ## ## $fitted.values ## [1] "numeric"
Hver værdi kan udtrækkes med $-tegnet efterfulgt af navnet på metrikken. For eksempel gemte du modellen som logit. For at udtrække AIC-kriterierne bruger du:
logit$aic
Output:
## [1] 27086.65
Trin 7) Vurder modellens ydeevne
Forvirringsmatrix
forvirringsmatrix er et bedre valg til at evaluere klassificeringsydelsen sammenlignet med de forskellige metrics, du så før. Den generelle idé er at tælle antallet af gange, hvor sande tilfælde klassificeres er falske.
For at beregne forvirringsmatricen skal du først have et sæt forudsigelser, så de kan sammenlignes med de faktiske mål.
predict <- predict(logit, data_test, type = 'response') # confusion matrix table_mat <- table(data_test$income, predict > 0.5) table_mat
Kode Forklaring
- predict(logit,data_test, type = 'respons'): Beregn forudsigelsen på testsættet. Sæt type = 'respons' for at beregne svarsandsynligheden.
- tabel(data_test$indkomst, forudsig > 0.5): Beregn forvirringsmatricen. forudsige > 0.5 betyder, at det returnerer 1, hvis de forudsagte sandsynligheder er over 0.5, ellers 0.
Output:
## ## FALSE TRUE ## <=50K 6310 495 ## >50K 1074 1229
Hver række i en forvirringsmatrix repræsenterer et faktisk mål, mens hver kolonne repræsenterer et forudsagt mål. Den første række af denne matrix betragter indkomsten lavere end 50k (den falske klasse): 6241 blev korrekt klassificeret som individer med en indkomst lavere end 50k (Sand negativ), mens den resterende fejlagtigt blev klassificeret som over 50k (Falsk positiv). Den anden række betragter indkomsten over 50k, den positive klasse var 1229 (Rigtig positiv), Mens Sand negativ var 1074.
Du kan beregne modellen nøjagtighed ved at summere det sande positive + det sande negative over den samlede observation
accuracy_Test <- sum(diag(table_mat)) / sum(table_mat) accuracy_Test
Kode Forklaring
- sum(diag(table_mat)): Summen af diagonalen
- sum(tabel_mat): Summen af matricen.
Output:
## [1] 0.8277339
Modellen ser ud til at lide af ét problem, den overvurderer antallet af falske negativer. Dette kaldes nøjagtighedstest paradoks. Vi sagde, at nøjagtigheden er forholdet mellem korrekte forudsigelser og det samlede antal tilfælde. Vi kan have relativt høj nøjagtighed, men en ubrugelig model. Det sker, når der er en dominerende klasse. Hvis du ser tilbage på forvirringsmatricen, kan du se, at de fleste tilfælde er klassificeret som sande negative. Forestil dig nu, at modellen klassificerede alle klasserne som negative (dvs. lavere end 50k). Du ville have en nøjagtighed på 75 procent (6718/6718+2257). Din model klarer sig bedre, men kæmper for at skelne det sande positive med det sande negative.
I en sådan situation er det at foretrække at have en mere kortfattet metrisk. Vi kan se på:
- Præcision=TP/(TP+FP)
- Recall=TP/(TP+FN)
Præcision vs Recall
Precision ser på nøjagtigheden af den positive forudsigelse. Recall er forholdet mellem positive tilfælde, der detekteres korrekt af klassificeringsorganet;
Du kan konstruere to funktioner til at beregne disse to metrikker
- Konstruer præcision
precision <- function(matrix) {
# True positive
tp <- matrix[2, 2]
# false positive
fp <- matrix[1, 2]
return (tp / (tp + fp))
}
Kode Forklaring
- mat[1,1]: Returner den første celle i den første kolonne i datarammen, dvs. den sande positive
- måtte[1,2]; Returner den første celle i den anden kolonne i datarammen, dvs. den falske positive
recall <- function(matrix) {
# true positive
tp <- matrix[2, 2]# false positive
fn <- matrix[2, 1]
return (tp / (tp + fn))
}
Kode Forklaring
- mat[1,1]: Returner den første celle i den første kolonne i datarammen, dvs. den sande positive
- måtte[2,1]; Returner den anden celle i den første kolonne i datarammen, dvs. det falske negativ
Du kan teste dine funktioner
prec <- precision(table_mat) prec rec <- recall(table_mat) rec
Output:
## [1] 0.712877 ## [2] 0.5336518
Når modellen siger, at det er et individ over 50k, er det kun korrekt i 54 procent af tilfældene og kan kræve personer over 50k i 72 procent af tilfældene.
Du kan oprette 
score baseret på præcision og genkaldelse. Det er et harmonisk gennemsnit af disse to metrikker, hvilket betyder, at det giver større vægt til de lavere værdier.
f1 <- 2 * ((prec * rec) / (prec + rec)) f1
Output:
## [1] 0.6103799
Afvejning mellem præcision og tilbagekaldelse
Det er umuligt at have både en høj præcision og høj genkaldelse.
Hvis vi øger præcisionen, vil det korrekte individ blive bedre forudsagt, men vi ville savne mange af dem (lavere tilbagekaldelse). I nogle situationer foretrækker vi højere præcision end tilbagekaldelse. Der er et konkavt forhold mellem præcision og genkaldelse.
- Forestil dig, du skal forudsige, om en patient har en sygdom. Du vil gerne være så præcis som muligt.
- Hvis du har brug for at opdage potentielle svigagtige personer på gaden gennem ansigtsgenkendelse, ville det være bedre at fange mange mennesker, der er stemplet som svigagtige, selvom præcisionen er lav. Politiet vil være i stand til at løslade den ikke-bedrageriske person.
ROC-kurven
Modtager Operating Karakteristisk kurve er et andet almindeligt værktøj, der bruges med binær klassificering. Den minder meget om præcision/genkaldelseskurven, men i stedet for at plotte præcision versus genkald, viser ROC-kurven den sande positive rate (dvs. genkaldelse) mod den falske positive rate. Den falske positive rate er forholdet mellem negative tilfælde, der er forkert klassificeret som positive. Det er lig med én minus den sande negative kurs. Den sande negative sats kaldes også specificitet. Derfor plotter ROC-kurven følsomhed (tilbagekaldelse) versus 1-specificitet
For at plotte ROC-kurven skal vi installere et bibliotek kaldet RORC. Vi kan finde i condaen bibliotek. Du kan indtaste koden:
conda installer -cr r-rocr -ja
Vi kan plotte ROC'en med funktionerne forudsigelse() og ydeevne().
library(ROCR) ROCRpred <- prediction(predict, data_test$income) ROCRperf <- performance(ROCRpred, 'tpr', 'fpr') plot(ROCRperf, colorize = TRUE, text.adj = c(-0.2, 1.7))
Kode Forklaring
- forudsigelse(forudsig, data_test$indkomst): ROCR-biblioteket skal oprette et forudsigelsesobjekt for at transformere inputdataene
- performance(ROCRpred, 'tpr','fpr'): Returner de to kombinationer, der skal produceres i grafen. Her er tpr og fpr konstrueret. For at plotte præcision og genkalde sammen, brug "prec", "rec".
Output:
Trin 8) Forbedre modellen
Du kan prøve at tilføje ikke-linearitet til modellen med interaktionen mellem
- alder og timer.om ugen
- køn og timer.om ugen.
Du skal bruge scoretesten til at sammenligne begge modeller
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
Output:
## [1] 0.6109181
Scoren er lidt højere end den forrige. Du kan fortsætte med at arbejde på dataene og prøve at slå resultatet.
Resumé
Vi kan opsummere funktionen til at træne en logistisk regression i nedenstående tabel:
| Pakke | Objektiv | Funktion | Argument |
|---|---|---|---|
| - | Opret tog-/testdatasæt | create_train_set() | data, størrelse, tog |
| GLM | Træn en generaliseret lineær model | glm() | formel, data, familie* |
| GLM | Opsummer modellen | Resumé() | monteret model |
| bund | Lav forudsigelse | forudsige() | tilpasset model, datasæt, type = 'svar' |
| bund | Lav en forvirringsmatrix | bord() | y, forudsige() |
| bund | Opret nøjagtighedsscore | sum(diag(tabel())/sum(tabel() | |
| ROCR | Opret ROC: Trin 1 Opret forudsigelse | forudsigelse() | forudsige(), y |
| ROCR | Opret ROC: Trin 2 Opret ydeevne | ydeevne() | forudsigelse(), 'tpr', 'fpr' |
| ROCR | Opret ROC: Trin 3 Plot graf | grund() | ydeevne() |
Den anden GLM type modeller er:
– binomial: (link = “logit”)
– gaussisk: (link = "identitet")
– Gamma: (link = "omvendt")
– inverse.gaussian: (link = “1/mu^2”)
– gift: (link = "log")
– kvasi: (link = "identitet", varians = "konstant")
– kvasibinomial: (link = "logit")
– quasipoisson: (link = "log")
