GLM in R: Verallgemeinertes lineares Modell mit Beispiel
Evelyn Clarke
Was ist logistische Regression?
Die logistische Regression wird verwendet, um eine Klasse, also eine Wahrscheinlichkeit, vorherzusagen. Die logistische Regression kann ein binรคres Ergebnis genau vorhersagen.
Stellen Sie sich vor, Sie mรถchten anhand vieler Merkmale vorhersagen, ob ein Kredit abgelehnt/angenommen wird. Die logistische Regression hat die Form 0/1. y = 0, wenn ein Kredit abgelehnt wird, y = 1, wenn er angenommen wird.
Ein logistisches Regressionsmodell unterscheidet sich vom linearen Regressionsmodell in zweierlei Hinsicht.
- Erstens akzeptiert die logistische Regression nur dichotomische (binรคre) Eingaben als abhรคngige Variable (dh einen Vektor von 0 und 1).
- Zweitens wird das Ergebnis durch die folgende probabilistische Linkfunktion gemessen, genannt Sigmoid aufgrund seiner S-Form.:
Die Ausgabe der Funktion liegt immer zwischen 0 und 1. Sehen Sie sich das Bild unten an
Die Sigmoidfunktion gibt Werte von 0 bis 1 zurรผck. Fรผr die Klassifizierungsaufgabe benรถtigen wir eine diskrete Ausgabe von 0 oder 1.
Um einen kontinuierlichen Fluss in einen diskreten Wert umzuwandeln, kรถnnen wir eine Entscheidungsgrenze von 0.5 festlegen. Alle Werte รผber diesem Schwellenwert werden als 1 klassifiziert
So erstellen Sie ein Generalized Liner Model (GLM)
Verwenden wir die Erwachsenen Datensatz zur Veranschaulichung der logistischen Regression. Der โErwachseneโ ist ein groรartiger Datensatz fรผr die Klassifizierungsaufgabe. Das Ziel besteht darin, vorherzusagen, ob das jรคhrliche Dollareinkommen einer Person 50.000 รผbersteigen wird. Der Datensatz enthรคlt 46,033 Beobachtungen und zehn Merkmale:
- Alter: Alter der Person. Numerisch
- Bildung: Bildungsniveau des Einzelnen. Faktor.
- Familienstand: MariGesamtstatus des Individuums. Faktor, z. B. nie verheiratet, verheirateter Lebenspartner, โฆ
- Geschlecht: Geschlecht der Person. Faktor, also mรคnnlich oder weiblich
- Einkommen: Target variabel. Einkommen รผber oder unter 50K. Faktor d. h. >50K, <=50K
unter anderem
library(dplyr)
data_adult <-read.csv("https://raw.githubusercontent.com/guru99-edu/R-Programming/master/adult.csv")
glimpse(data_adult)
Ausgang:
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...
Wir werden wie folgt vorgehen:
- Schritt 1: Kontinuierliche Variablen prรผfen
- Schritt 2: Faktorvariablen prรผfen
- Schritt 3: Feature-Engineering
- Schritt 4: Zusammenfassende Statistik
- Schritt 5: Trainings-/Testset
- Schritt 6: Erstellen Sie das Modell
- Schritt 7: Bewerten Sie die Leistung des Modells
- Schritt 8: Verbessern Sie das Modell
Ihre Aufgabe besteht darin, vorherzusagen, welche Person einen Umsatz von mehr als 50 erzielen wird.
In diesem Tutorial wird jeder Schritt detailliert beschrieben, um eine Analyse an einem realen Datensatz durchzufรผhren.
Schritt 1) โโรberprรผfen Sie kontinuierliche Variablen
Im ersten Schritt sehen Sie die Verteilung der kontinuierlichen Variablen.
continuous <-select_if(data_adult, is.numeric) summary(continuous)
Code Erlรคuterung
- kontinuierlich <- select_if(data_adult, is.numeric): Verwenden Sie die Funktion select_if() aus der dplyr-Bibliothek, um nur die numerischen Spalten auszuwรคhlen
- Zusammenfassung (kontinuierlich): Drucken Sie die zusammenfassende Statistik
Ausgang:
## 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
Aus der obigen Tabelle kรถnnen Sie ersehen, dass die Daten vรถllig unterschiedliche Skalen aufweisen und dass es bei den Stunden pro Woche groรe Ausreiรer gibt (sehen Sie sich also das letzte Quartil und den Maximalwert an).
Sie kรถnnen das Problem in den folgenden zwei Schritten lรถsen:
- 1: Zeichnen Sie die Verteilung der Stunden pro Woche auf
- 2: Standardisieren Sie die kontinuierlichen Variablen
- Plotten Sie die Verteilung
Schauen wir uns die Verteilung der Stunden pro Woche genauer an
# Histogram with kernel density curve
library(ggplot2)
ggplot(continuous, aes(x = hours.per.week)) +
geom_density(alpha = .2, fill = "#FF6666")
Ausgang:
Die Variable hat viele Ausreiรer und eine nicht gut definierte Verteilung. Sie kรถnnen dieses Problem teilweise lรถsen, indem Sie die oberen 0.01 Prozent der Stunden pro Woche lรถschen.
Grundlegende Syntax des Quantils:
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.
Wir berechnen das oberste 2-Prozent-Perzentil
top_one_percent <- quantile(data_adult$hours.per.week, .99) top_one_percent
Code Erlรคuterung
- quantile(data_adult$hours.per.week, .99): Berechnen Sie den Wert der 99 Prozent der Arbeitszeit
Ausgang:
## 99% ## 80
98 Prozent der Bevรถlkerung arbeiten weniger als 80 Stunden pro Woche.
Sie kรถnnen die Beobachtungen รผber diesem Schwellenwert lรถschen. Sie verwenden den Filter aus dem dplyr Bibliothek.
data_adult_drop <-data_adult %>% filter(hours.per.week<top_one_percent) dim(data_adult_drop)
Ausgang:
## [1] 45537 10
- Standardisieren Sie die kontinuierlichen Variablen
Sie kรถnnen jede Spalte standardisieren, um die Leistung zu verbessern, da Ihre Daten nicht den gleichen Maรstab haben. Sie kรถnnen die Funktion mutate_if aus der dplyr-Bibliothek verwenden. Die grundlegende Syntax lautet:
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
Sie kรถnnen die numerischen Spalten wie folgt standardisieren:
data_adult_rescale <- data_adult_drop % > % mutate_if(is.numeric, funs(as.numeric(scale(.)))) head(data_adult_rescale)
Code Erlรคuterung
- mutate_if(is.numeric, funs(scale)): Die Bedingung ist nur eine numerische Spalte und die Funktion ist skaliert
Ausgang:
## 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
Schritt 2) รberprรผfen Sie die Faktorvariablen
Dieser Schritt hat zwei Ziele:
- รberprรผfen Sie die Ebene in jeder kategorialen Spalte
- Definieren Sie neue Ebenen
Wir werden diesen Schritt in drei Teile unterteilen:
- Wรคhlen Sie die kategorialen Spalten aus
- Speichern Sie das Balkendiagramm jeder Spalte in einer Liste
- Drucken Sie die Grafiken aus
Wir kรถnnen die Faktorspalten mit dem folgenden Code auswรคhlen:
# Select categorical column factor <- data.frame(select_if(data_adult_rescale, is.factor)) ncol(factor)
Code Erlรคuterung
- data.frame(select_if(data_adult, is.factor)): Wir speichern die Faktorspalten in Faktor in einem Datenrahmentyp. Die Bibliothek ggplot2 erfordert ein Datenrahmenobjekt.
Ausgang:
## [1] 6
Der Datensatz enthรคlt 6 kategoriale Variablen
Der zweite Schritt ist geschickter. Sie mรถchten fรผr jede Spalte im Datenrahmenfaktor ein Balkendiagramm zeichnen. Es ist bequemer, den Prozess zu automatisieren, insbesondere wenn viele Spalten vorhanden sind.
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)))
Code Erlรคuterung
- lapply(): Verwenden Sie die Funktion lapply(), um eine Funktion in allen Spalten des Datensatzes zu รผbergeben. Sie speichern die Ausgabe in einer Liste
- Funktion(x): Die Funktion wird fรผr jedes x abgearbeitet. Hier sind x die Spalten
- ggplot(factor, aes(get(x))) + geom_bar()+ theme(axis.text.x = element_text(angle = 90)): Erstellen Sie ein Balkendiagramm fรผr jedes x-Element. Hinweis: Um x als Spalte zurรผckzugeben, mรผssen Sie es in get() einschlieรen.
Der letzte Schritt ist relativ einfach. Sie mรถchten die 6 Grafiken ausdrucken.
# Print the graph graph
Ausgang:
## [[1]]
## ## [[2]]
## ## [[3]]
## ## [[4]]
## ## [[5]]
## ## [[6]]
Hinweis: Verwenden Sie die Schaltflรคche โWeiterโ, um zum nรคchsten Diagramm zu navigieren
Schritt 3) Feature-Engineering
Neufassung der Bildung
Aus der obigen Grafik kรถnnen Sie ersehen, dass die Variable Bildung 16 Stufen hat. Dies ist erheblich, und einige Ebenen weisen eine relativ geringe Anzahl von Beobachtungen auf. Wenn Sie die Menge an Informationen, die Sie aus dieser Variable erhalten kรถnnen, verbessern mรถchten, kรถnnen Sie sie auf eine hรถhere Ebene umwandeln. Sie bilden nรคmlich grรถรere Gruppen mit รคhnlichem Bildungsniveau. Beispielsweise wird ein niedriges Bildungsniveau in einen Schulabbruch umgewandelt. Hรถhere Bildungsstufen werden auf Master umgestellt.
Hier ist das Detail:
| Altes Niveau | Neues level |
|---|---|
| vorschulisch | aussteigen |
| 10. | Aussteiger |
| 11. | Aussteiger |
| 12. | Aussteiger |
| 1.-4 | Aussteiger |
| 5th-6th | Aussteiger |
| 7th-8th | Aussteiger |
| 9. | Aussteiger |
| HS-Grad | HighGrad |
| Einige College | Gemeinschaft |
| Assoc-acdm | Gemeinschaft |
| Assoc-voc | Gemeinschaft |
| Junggesellen | Junggesellen |
| Masters | Masters |
| Prof-Schule | Masters |
| Doktor | 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")))))))
Code Erlรคuterung
- Wir verwenden das Verb mutate aus der dplyr-Bibliothek. Mit der Aussage ifelse verรคndern wir die Werte der Bildung
In der folgenden Tabelle erstellen Sie eine zusammenfassende Statistik, um zu sehen, wie viele Ausbildungsjahre (Z-Wert) durchschnittlich erforderlich sind, um den Bachelor, Master oder Ph.D. zu erreichen.
recast_data % > % group_by(education) % > % summarize(average_educ_year = mean(educational.num), count = n()) % > % arrange(average_educ_year)
Ausgang:
## # 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
Neufassung MariTal-Status
Es ist auch mรถglich, niedrigere Ebenen fรผr den Familienstand zu erstellen. Im folgenden Code รคndern Sie die Ebene wie folgt:
| Altes Niveau | Neues level |
|---|---|
| Nie verheiratet | Nicht verheiratet |
| Verheirateter-Ehepartner-abwesend | Nicht verheiratet |
| Verheirateter AF-Ehepartner | Verheiratet |
| Verheirateter-Lebenspartner | |
| Getrennt | Getrennt |
| Geschieden | |
| Witwen | Witwe |
# 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")))))
Sie kรถnnen die Anzahl der Personen innerhalb jeder Gruppe รผberprรผfen.
table(recast_data$marital.status)
Ausgang:
## ## Married Not_married Separated Widow ## 21165 15359 7727 1286
Schritt 4) Zusammenfassende Statistik
Es ist an der Zeit, einige Statistiken zu unseren Zielvariablen zu รผberprรผfen. In der Grafik unten zรคhlen Sie den Prozentsatz der Personen, die aufgrund ihres Geschlechts mehr als 50 verdienen.
# Plot gender income
ggplot(recast_data, aes(x = gender, fill = income)) +
geom_bar(position = "fill") +
theme_classic()
Ausgang:
Prรผfen Sie als Nรคchstes, ob die Herkunft der Person Einfluss auf ihr Einkommen hat.
# 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))
Ausgang:
Die Anzahl der Arbeitszeiten nach Geschlecht.
# 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()
Ausgang:
Der Boxplot bestรคtigt, dass die Arbeitszeitverteilung auf unterschiedliche Gruppen passt. Im Boxplot gibt es keine homogenen Beobachtungen beider Geschlechter.
Sie kรถnnen die Dichte der wรถchentlichen Arbeitszeit nach Bildungsart รผberprรผfen. Die Verteilungen weisen viele unterschiedliche Ausprรคgungen auf. Dies lรคsst sich wahrscheinlich durch die Art der Ausbildung erklรคren.tract in den 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()
Code Erlรคuterung
- ggplot(recast_data, aes( x= hours.per.week)): Ein Dichtediagramm erfordert nur eine Variable
- geom_density(aes(color = education), alpha =0.5): Das geometrische Objekt zur Steuerung der Dichte
Ausgang:
Um Ihre Gedanken zu bestรคtigen, kรถnnen Sie eine Einbahnstraรe durchfรผhren ANOVA-Test:
anova <- aov(hours.per.week~education, recast_data) summary(anova)
Ausgang:
## 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
Der ANOVA-Test bestรคtigt den Unterschied im Durchschnitt zwischen den Gruppen.
Nichtlinearitรคt
Bevor Sie das Modell ausfรผhren, kรถnnen Sie prรผfen, ob die Anzahl der geleisteten Arbeitsstunden mit dem Alter zusammenhรคngt.
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()
Code Erlรคuterung
- ggplot(recast_data, aes(x = age, y = hours.per.week)): Legt die รsthetik des Diagramms fest
- geom_point(aes(color= Einkommen), Grรถรe =0.5): Konstruieren Sie das Punktdiagramm
- stat_smooth(): Fรผgen Sie die Trendlinie mit den folgenden Argumenten hinzu:
- method='lm': Zeichnen Sie den angepassten Wert auf, wenn der lineare Regression
- Formel = y~poly(x,2): Passen Sie eine polynomielle Regression an
- se = TRUE: Standardfehler hinzufรผgen
- aes(color= Einkommen): Brechen Sie das Modell nach Einkommen auf
Ausgang:
Kurz gesagt: Sie kรถnnen Interaktionsterme im Modell testen, um den Nichtlinearitรคtseffekt zwischen der wรถchentlichen Arbeitszeit und anderen Merkmalen zu erfassen. Es ist wichtig zu erkennen, unter welchen Bedingungen sich die Arbeitszeit unterscheidet.
Korrelation
Die nรคchste Prรผfung besteht darin, die Korrelation zwischen den Variablen zu visualisieren. Sie wandeln den Typ der Faktorebene in einen numerischen Typ um, sodass Sie eine Wรคrmekarte zeichnen kรถnnen, die den mit der Spearman-Methode berechneten Korrelationskoeffizienten enthรคlt.
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")
Code Erlรคuterung
- data.frame(lapply(recast_data,as.integer)): Daten in numerisch konvertieren
- ggcorr() zeichnet die Heatmap mit den folgenden Argumenten:
- Methode: Methode zur Berechnung der Korrelation
- nbreaks = 6: Anzahl der Pausen
- hjust = 0.8: Kontrollposition des Variablennamens im Plot
- label = TRUE: Beschriftungen in der Mitte der Fenster hinzufรผgen
- label_size = 3: Grรถรenbeschriftungen
- color = โgrey50โ): Farbe des Etiketts
Ausgang:
Schritt 5) Trainings-/Testset
Jeder wird beaufsichtigt Maschinelles Lernen Fรผr diese Aufgabe ist es erforderlich, die Daten zwischen einem Zugsatz und einem Testsatz aufzuteilen. Sie kรถnnen die โFunktionโ, die Sie in den anderen Tutorials zum betreuten Lernen erstellt haben, verwenden, um einen Zug-/Testsatz zu erstellen.
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)
Ausgang:
## [1] 36429 9
dim(data_test)
Ausgang:
## [1] 9108 9
Schritt 6) Erstellen Sie das Modell
Um zu sehen, wie der Algorithmus funktioniert, verwenden Sie das Paket glm(). Der Verallgemeinertes lineares Modell ist eine Sammlung von Modellen. Die grundlegende Syntax lautet:
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")
Sie sind bereit, das Logistikmodell zu schรคtzen, um das Einkommensniveau auf eine Reihe von Merkmalen aufzuteilen.
formula <- income~. logit <- glm(formula, data = data_train, family = 'binomial') summary(logit)
Code Erlรคuterung
- Formel <- Einkommen ~ .: Erstellen Sie das passende Modell
- logit <- glm(formula, data = data_train, family = 'binomial'): Passen Sie ein Logistikmodell (family = 'binomial') an die data_train-Daten an.
- summary(logit): Drucken Sie die Zusammenfassung des Modells
Ausgang:
## ## 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
Die Zusammenfassung unseres Modells enthรผllt interessante Informationen. Die Leistung einer logistischen Regression wird anhand spezifischer Schlรผsselmetriken bewertet.
- AIC (Akaike Information Criteria): Dies ist das รquivalent von R2 in der logistischen Regression. Es misst die Anpassung, wenn eine Strafe auf die Anzahl der Parameter angewendet wird. Kleiner AIC Werte zeigen an, dass das Modell nรคher an der Wahrheit liegt.
- Nullabweichung: Passt das Modell nur mit dem Achsenabschnitt an. Der Freiheitsgrad betrรคgt n-1. Wir kรถnnen es als Chi-Quadrat-Wert interpretieren (angepasster Wert, der sich vom tatsรคchlichen Wert der Hypothesenprรผfung unterscheidet).
- Restabweichung: Modell mit allen Variablen. Es wird auch als Chi-Quadrat-Hypothesetest interpretiert.
- Anzahl der Fisher-Scoring-Iterationen: Anzahl der Iterationen vor der Konvergenz.
Die Ausgabe der Funktion glm() wird in einer Liste gespeichert. Der folgende Code zeigt alle Elemente, die in der Logit-Variablen verfรผgbar sind, die wir zur Auswertung der logistischen Regression erstellt haben.
# Die Liste ist sehr lang. Geben Sie nur die ersten drei Elemente aus
lapply(logit, class)[1:3]
Ausgang:
## $coefficients ## [1] "numeric" ## ## $residuals ## [1] "numeric" ## ## $fitted.values ## [1] "numeric"
Jeder Wert kann ex seintracted mit dem $-Zeichen gefolgt vom Namen der Metriken. Zum Beispiel haben Sie das Modell als logit gespeichert. Um z. B.tracBei den AIC-Kriterien verwenden Sie:
logit$aic
Ausgang:
## [1] 27086.65
Schritt 7) Bewerten Sie die Leistung des Modells
Verwirrung Matrix
Das Verwirrung Matrix ist eine bessere Wahl zur Bewertung der Klassifizierungsleistung im Vergleich zu den verschiedenen Metriken, die Sie zuvor gesehen haben. Die allgemeine Idee besteht darin, zu zรคhlen, wie oft Wahre Instanzen als Falsch klassifiziert werden.
Um die Verwirrungsmatrix zu berechnen, benรถtigen Sie zunรคchst eine Reihe von Vorhersagen, damit diese mit den tatsรคchlichen Zielen verglichen werden kรถnnen.
predict <- predict(logit, data_test, type = 'response') # confusion matrix table_mat <- table(data_test$income, predict > 0.5) table_mat
Code Erlรคuterung
- Predict(logit,data_test, type = 'response'): Berechnen Sie die Vorhersage fรผr den Testsatz. Legen Sie type = 'response' fest, um die Antwortwahrscheinlichkeit zu berechnen.
- table(data_test$income, Predict > 0.5): Berechnen Sie die Verwirrungsmatrix. Vorhersage > 0.5 bedeutet, dass 1 zurรผckgegeben wird, wenn die vorhergesagten Wahrscheinlichkeiten รผber 0.5 liegen, andernfalls 0.
Ausgang:
## ## FALSE TRUE ## <=50K 6310 495 ## >50K 1074 1229
Jede Zeile in einer Verwirrungsmatrix stellt ein tatsรคchliches Ziel dar, wรคhrend jede Spalte ein vorhergesagtes Ziel darstellt. Die erste Zeile dieser Matrix berรผcksichtigt das Einkommen von weniger als 50 (die Falsch-Klasse): 6241 wurden korrekt als Personen mit einem Einkommen von weniger als 50 klassifiziert (Wahr-negativ), wรคhrend der verbleibende fรคlschlicherweise als รผber 50 eingestuft wurde (Falsch positiv). Die zweite Zeile berรผcksichtigt das Einkommen รผber 50, die positive Klasse war 1229 (Richtig positiv), Wรคhrend die Wahr-negativ war 1074.
Sie kรถnnen das Modell berechnen Genauigkeit durch Summieren des wahren Positivs + des wahren Negativs รผber die Gesamtbeobachtung
accuracy_Test <- sum(diag(table_mat)) / sum(table_mat) accuracy_Test
Code Erlรคuterung
- sum(diag(table_mat)): Summe der Diagonalen
- sum(table_mat): Summe der Matrix.
Ausgang:
## [1] 0.8277339
Das Modell scheint unter einem Problem zu leiden: Es รผberschรคtzt die Anzahl der falsch-negativen Ergebnisse. Dies nennt man Genauigkeitstest-Paradoxon. Wir haben angegeben, dass die Genauigkeit das Verhรคltnis der richtigen Vorhersagen zur Gesamtzahl der Fรคlle ist. Wir kรถnnen eine relativ hohe Genauigkeit, aber ein nutzloses Modell haben. Es passiert, wenn es eine dominierende Klasse gibt. Wenn Sie sich die Verwirrungsmatrix noch einmal ansehen, kรถnnen Sie sehen, dass die meisten Fรคlle als richtig negativ eingestuft werden. Stellen Sie sich nun vor, das Modell hรคtte alle Klassen als negativ klassifiziert (dh niedriger als 50). Sie hรคtten eine Genauigkeit von 75 Prozent (6718/6718+2257). Ihr Modell bietet eine bessere Leistung, hat jedoch Schwierigkeiten, das wahre Positive vom wahren Negativen zu unterscheiden.
In einer solchen Situation ist es vorzuziehen, eine prรคzisere Metrik zu verwenden. Wir kรถnnen uns Folgendes ansehen:
- Prรคzision=TP/(TP+FP)
- Rรผckruf=TP/(TP+FN)
Prรคzision vs. Rรผckruf
Prรคzision untersucht die Genauigkeit der positiven Vorhersage. Erinnern ist das Verhรคltnis der positiven Instanzen, die vom Klassifikator korrekt erkannt werden;
Sie kรถnnen zwei Funktionen erstellen, um diese beiden Metriken zu berechnen
- Konstruieren Sie Prรคzision
precision <- function(matrix) {
# True positive
tp <- matrix[2, 2]
# false positive
fp <- matrix[1, 2]
return (tp / (tp + fp))
}
Code Erlรคuterung
- mat[1,1]: Gibt die erste Zelle der ersten Spalte des Datenrahmens zurรผck, also das wahre Positive
- mat[1,2]; Gibt die erste Zelle der zweiten Spalte des Datenrahmens zurรผck, also das falsche Positiv
recall <- function(matrix) {
# true positive
tp <- matrix[2, 2]# false positive
fn <- matrix[2, 1]
return (tp / (tp + fn))
}
Code Erlรคuterung
- mat[1,1]: Gibt die erste Zelle der ersten Spalte des Datenrahmens zurรผck, also das wahre Positive
- mat[2,1]; Gibt die zweite Zelle der ersten Spalte des Datenrahmens zurรผck, also das falsche Negativ
Sie kรถnnen Ihre Funktionen testen
prec <- precision(table_mat) prec rec <- recall(table_mat) rec
Ausgang:
## [1] 0.712877 ## [2] 0.5336518
Wenn das Modell angibt, dass es sich um eine Person รผber 50 handelt, ist es nur in 54 Prozent der Fรคlle korrekt und kann in 50 Prozent der Fรคlle Personen รผber 72 angeben.
Sie kรถnnen die erstellen 
Bewertung basierend auf Prรคzision und Erinnerung. Der ist ein harmonisches Mittel dieser beiden Metriken, was bedeutet, dass den niedrigeren Werten mehr Gewicht beigemessen wird.
f1 <- 2 * ((prec * rec) / (prec + rec)) f1
Ausgang:
## [1] 0.6103799
Kompromiss zwischen Prรคzision und Rรผckruf
Es ist unmรถglich, gleichzeitig eine hohe Prรคzision und einen hohen Rรผckruf zu erreichen.
Wenn wir die Prรคzision erhรถhen, wird die richtige Person besser vorhergesagt, wir wรผrden jedoch viele davon รผbersehen (geringere Erinnerung). In manchen Situationen bevorzugen wir eine hรถhere Prรคzision als einen Rรผckruf. Es besteht ein konkaver Zusammenhang zwischen Prรคzision und Erinnerung.
- Stellen Sie sich vor, Sie mรผssen vorhersagen, ob ein Patient eine Krankheit hat. Sie mรถchten so prรคzise wie mรถglich sein.
- Wenn Sie potenziell betrรผgerische Personen auf der Straรe mithilfe der Gesichtserkennung erkennen mรผssen, wรคre es besser, viele als betrรผgerisch eingestufte Personen zu erfassen, auch wenn die Prรคzision gering ist. Die Polizei kann die nicht betrรผgerische Person freilassen.
Die ROC-Kurve
Das Empfรคnger Operating-Charakteristik Kurve ist ein weiteres gรคngiges Werkzeug, das bei der binรคren Klassifizierung verwendet wird. Sie ist der Prรคzisions-/Erinnerungskurve sehr รคhnlich, aber anstatt Prรคzision gegen Erinnerung darzustellen, zeigt die ROC-Kurve die wahre positive Rate (d. h. Erinnerung) gegenรผber der falsch positiven Rate. Die Falsch-Positiv-Rate ist das Verhรคltnis der negativen Fรคlle, die fรคlschlicherweise als positiv eingestuft werden. Er entspricht eins minus der echten Negativrate. Man spricht auch von der True-Negativ-Rate Spezifitรคt. Daher die ROC-Kurvendiagramme Empfindlichkeit (Rรผckruf) versus 1-Spezifitรคt
Um die ROC-Kurve darzustellen, mรผssen wir eine Bibliothek namens RORC installieren. Wir kรถnnen in der Conda finden . Sie kรถnnen den Code eingeben:
conda install -cr r-rocr โyes
Wir kรถnnen den ROC mit den Funktionen โprediction()โ und โperformance()โ darstellen.
library(ROCR) ROCRpred <- prediction(predict, data_test$income) ROCRperf <- performance(ROCRpred, 'tpr', 'fpr') plot(ROCRperf, colorize = TRUE, text.adj = c(-0.2, 1.7))
Code Erlรคuterung
- Vorhersage (Vorhersage, Datentest $ Einkommen): Die ROCR-Bibliothek muss ein Vorhersageobjekt erstellen, um die Eingabedaten zu transformieren
- performance(ROCRpred, 'tpr','fpr'): Gibt die beiden Kombinationen zurรผck, die im Diagramm erzeugt werden sollen. Hier werden tpr und fpr konstruiert. Um die Plotgenauigkeit und den Abruf zusammen zu erreichen, verwenden Sie โprecโ und โrecโ.
Ausgang:
Schritt 8) Verbessern Sie das Modell
Sie kรถnnen versuchen, dem Modell durch die Interaktion zwischen Nichtlinearitรคt hinzuzufรผgen
- Alter und Stunden pro Woche
- Geschlecht und Stunden pro Woche.
Sie mรผssen den Score-Test verwenden, um beide Modelle zu vergleichen
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
Ausgang:
## [1] 0.6109181
Die Punktzahl ist etwas hรถher als die vorherige. Sie kรถnnen weiter an den Daten arbeiten und versuchen, die Punktzahl zu รผbertreffen.
Zusammenfassung
Wir kรถnnen die Funktion zum Trainieren einer logistischen Regression in der folgenden Tabelle zusammenfassen:
| Verpackung | Ziel | Funktion | Argument |
|---|---|---|---|
| - | Erstellen Sie einen Trainings-/Testdatensatz | create_train_set() | Daten, Grรถรe, Zug |
| glm | Trainieren Sie ein verallgemeinertes lineares Modell | glm() | Formel, Daten, Familie* |
| glm | Das Modell zusammenfassen | Zusammenfassung() | angepasstes Modell |
| Base | Machen Sie eine Vorhersage | vorhersagen() | angepasstes Modell, Datensatz, Typ = 'Antwort' |
| Base | Erstellen Sie eine Verwirrungsmatrix | Tisch() | y, vorhersagen() |
| Base | Genauigkeitswert erstellen | sum(diag(table())/sum(table() | |
| ROCR | ROC erstellen: Schritt 1 Vorhersage erstellen | Vorhersage() | vorhersagen(), y |
| ROCR | ROC erstellen: Schritt 2 Leistung erstellen | Leistung() | Vorhersage(), 'tpr', 'fpr' |
| ROCR | ROC erstellen: Schritt 3 Diagramm zeichnen | Handlung() | Leistung() |
Die andere GLM Art der Modelle sind:
โ Binomial: (link = โlogitโ)
โ Gauรsche Funktion: (Link = โIdentitรคtโ)
โ Gamma: (Link = โinversโ)
โ inverse.gaussian: (link = โ1/mu^2โ)
โ Poisson: (Link = โlogโ)
โ quasi: (Link = โIdentitรคtโ, Varianz = โKonstanteโ)
โ quasibinomial: (link = โlogitโ)
โ quasipoisson: (link = โlogโ)
