GLM in R: Mô hình tuyến tính tổng quát với ví dụ

Hồi quy Logistic là gì?

Hồi quy logistic được sử dụng để dự đoán một lớp, tức là một xác suất. Hồi quy logistic có thể dự đoán chính xác kết quả nhị phân.

Hãy tưởng tượng bạn muốn dự đoán liệu một khoản vay có bị từ chối / chấp nhận hay không dựa trên nhiều thuộc tính. Hồi quy logistic có dạng 0/1. y = 0 nếu khoản vay bị từ chối, y = 1 nếu được chấp nhận.

Mô hình hồi quy logistic khác với mô hình hồi quy tuyến tính theo hai cách.

• Trước hết, hồi quy logistic chỉ chấp nhận đầu vào lưỡng phân (nhị phân) như một biến phụ thuộc (tức là vectơ 0 và 1).
• Thứ hai, kết quả được đo bằng hàm liên kết xác suất sau đây được gọi là sigmoid do có hình chữ S.

Đầu ra của hàm luôn nằm trong khoảng từ 0 đến 1. Kiểm tra Hình ảnh bên dưới

Hàm sigmoid trả về các giá trị từ 0 đến 1. Đối với nhiệm vụ phân loại, chúng ta cần một đầu ra rời rạc là 0 hoặc 1.

Để chuyển đổi một luồng liên tục thành giá trị rời rạc, chúng ta có thể đặt một quyết định bị ràng buộc ở mức 0,5. Tất cả các giá trị trên ngưỡng này được phân loại là 1

Trong hướng dẫn này, bạn sẽ học

• Hồi quy Logistic là gì?
• Cách tạo Mô hình lót tổng quát (GLM)
• Bước 1) Kiểm tra các biến liên tục
• Bước 2) Kiểm tra các biến nhân tố
• Bước 3) Kỹ thuật tính năng
• Bước 4) Thống kê Tóm tắt
• Bước 5) Đào tạo / bộ kiểm tra
• Bước 6) Xây dựng mô hình
• Bước 7) Đánh giá hiệu suất của mô hình

Cách tạo Mô hình lót tổng quát (GLM)

Hãy sử dụng tập dữ liệu người lớn để minh họa hồi quy Logistic. "Người lớn" là một tập dữ liệu tuyệt vời cho nhiệm vụ phân loại. Mục tiêu là để dự đoán liệu thu nhập hàng năm tính bằng đô la của một cá nhân có vượt quá 50.000 hay không. Tập dữ liệu chứa 46.033 quan sát và mười tính năng:

• age: tuổi của cá nhân. Số
• học vấn: Trình độ học vấn của cá nhân. Hệ số.
• marital.status: Tình trạng hôn nhân của cá nhân. Yếu tố tức là Chưa từng kết hôn, Đã kết hôn-công dân-vợ / chồng,…
• giới tính: Giới tính của cá nhân. Yếu tố, tức là Nam hoặc Nữ
• thu nhập: Biến mục tiêu. Thu nhập trên hoặc dưới 50K. Hệ số tức là> 50K, <= 50K

giữa những người khác

library(dplyr)data_adult <-read.csv("https://raw.githubusercontent.com/guru99-edu/R-Programming/master/adult.csv")glimpse(data_adult)

Đầu ra:

Observations: 48,842Variables: 10$ x  1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,… $ age  25, 38, 28, 44, 18, 34, 29, 63, 24, 55, 65, 36, 26… $ workclass  Private, Private, Local-gov, Private, ?, Private,… $ education  11th, HS-grad, Assoc-acdm, Some-college, Some-col… $ educational.num  7, 9, 12, 10, 10, 6, 9, 15, 10, 4, 9, 13, 9, 9, 9,… $ marital.status  Never-married, Married-civ-spouse, Married-civ-sp… $ race  Black, White, White, Black, White, White, Black,… $ gender  Male, Male, Male, Male, Female, Male, Male, Male,… $ hours.per.week  40, 50, 40, 40, 30, 30, 40, 32, 40, 10, 40, 40, 39… $ income  <=50K, <=50K, >50K, >50K, <=50K, <=50K, <=50K, >5… 

Chúng tôi sẽ tiến hành như sau:

• Bước 1: Kiểm tra các biến liên tục
• Bước 2: Kiểm tra các biến nhân tố
• Bước 3: Kỹ thuật tính năng
• Bước 4: Thống kê tóm tắt
• Bước 5: Đào tạo / tập kiểm tra
• Bước 6: Xây dựng mô hình
• Bước 7: Đánh giá hiệu suất của mô hình
• bước 8: Cải thiện mô hình

Nhiệm vụ của bạn là dự đoán cá nhân nào sẽ có doanh thu cao hơn 50K.

Trong hướng dẫn này, từng bước sẽ được trình bày chi tiết để thực hiện phân tích trên tập dữ liệu thực.

Bước 1) Kiểm tra các biến liên tục

Trong bước đầu tiên, bạn có thể thấy sự phân bố của các biến liên tục.

continuous <-select_if(data_adult, is.numeric)summary(continuous)

Giải thích mã

• liên tục <- select_if (data_adult, is.numeric): Sử dụng hàm select_if () từ thư viện dplyr để chỉ chọn các cột số
• tóm tắt (liên tục): In thống kê tóm tắt

Đầu ra:

## 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

Từ bảng trên, bạn có thể thấy rằng dữ liệu có quy mô và giờ hoàn toàn khác nhau .per.weeks có giá trị ngoại lệ lớn (.ie hãy nhìn vào phần tư cuối cùng và giá trị tối đa).

Bạn có thể giải quyết nó theo hai bước:

• 1: Vẽ biểu đồ phân phối giờ.per.week
• 2: Chuẩn hóa các biến liên tục

1. Lập đồ thị phân phối

Chúng ta hãy xem xét kỹ hơn sự phân bổ của giờ.per.week

# Histogram with kernel density curvelibrary(ggplot2)ggplot(continuous, aes(x = hours.per.week)) +geom_density(alpha = .2, fill = "#FF6666")

Đầu ra:

Biến có rất nhiều giá trị ngoại lai và phân phối không được xác định rõ. Bạn có thể giải quyết một phần vấn đề này bằng cách xóa 0,01% số giờ hàng tuần.

Cú pháp cơ bản của tập tin lượng tử:

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.

Chúng tôi tính toán 2 phần trăm hàng đầu

top_one_percent <- quantile(data_adult$hours.per.week, .99)top_one_percent

Giải thích mã

• quantile (data_adult $ hours.per.week, .99): Tính giá trị của 99 phần trăm thời gian làm việc

Đầu ra:

## 99%## 80 

98 phần trăm dân số làm việc dưới 80 giờ mỗi tuần.

Bạn có thể giảm các quan sát trên ngưỡng này. Bạn sử dụng bộ lọc từ thư viện dplyr.

data_adult_drop <-data_adult %>%filter(hours.per.week
Đầu ra:


## [1] 45537 10 
 

Chuẩn hóa các biến liên tục

Bạn có thể chuẩn hóa từng cột để cải thiện hiệu suất vì dữ liệu của bạn không có cùng tỷ lệ. Bạn có thể sử dụng hàm mutate_if từ thư viện dplyr. Cú pháp cơ bản là:
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
Bạn có thể chuẩn hóa các cột số như sau:
data_adult_rescale <- data_adult_drop % > %mutate_if(is.numeric, funs(as.numeric(scale(.))))head(data_adult_rescale)
Giải thích mã

mutate_if (is.numeric, funs (scale)): Điều kiện chỉ là cột số và hàm là tỷ lệ

Đầu ra:
## 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
Bước 2) Kiểm tra các biến nhân tố
Bước này có hai mục tiêu:

Kiểm tra mức độ trong mỗi cột phân loại
Xác định cấp độ mới

Chúng tôi sẽ chia bước này thành ba phần:

Chọn các cột phân loại
Lưu trữ biểu đồ thanh của mỗi cột trong danh sách
In các đồ thị

Chúng ta có thể chọn các cột yếu tố với mã bên dưới:
# Select categorical columnfactor <- data.frame(select_if(data_adult_rescale, is.factor))ncol(factor)
Giải thích mã

data.frame (select_if (data_adult, is.factor)): Chúng tôi lưu trữ các cột yếu tố trong yếu tố trong một kiểu khung dữ liệu. Thư viện ggplot2 yêu cầu một đối tượng khung dữ liệu.

Đầu ra:
## [1] 6 
Tập dữ liệu chứa 6 biến phân loại
Bước thứ hai là có kỹ năng hơn. Bạn muốn vẽ biểu đồ thanh cho mỗi cột trong hệ số khung dữ liệu. Sẽ thuận tiện hơn khi tự động hóa quy trình, đặc biệt là trong tình huống có rất nhiều cột.
library(ggplot2)# Create graph for each columngraph <- lapply(names(factor),function(x)ggplot(factor, aes(get(x))) +geom_bar() +theme(axis.text.x = element_text(angle = 90)))
Giải thích mã

lapply (): Sử dụng hàm lapply () để chuyển một hàm trong tất cả các cột của tập dữ liệu. Bạn lưu trữ kết quả đầu ra trong một danh sách
function (x): Hàm sẽ được xử lý cho mỗi x. Đây là cột
ggplot (factor, aes (get (x))) + geom_bar () + theme (axis.text.x = element_text (angle = 90)): Tạo biểu đồ thanh char cho mỗi phần tử x. Lưu ý, để trả về x dưới dạng một cột, bạn cần đưa nó vào bên trong get ()

Bước cuối cùng tương đối dễ dàng. Bạn muốn in 6 đồ thị.
# Print the graphgraph
Đầu ra:
## [[1]]

## ## [[2]]

## ## [[3]]

## ## [[4]]

## ## [[5]]

## ## [[6]]

Lưu ý: Sử dụng nút tiếp theo để điều hướng đến biểu đồ tiếp theo

Bước 3) Kỹ thuật tính năng
Đúc kết lại giáo dục
Từ biểu đồ trên, bạn có thể thấy rằng trình độ học vấn có 16 cấp độ. Điều này là đáng kể và một số cấp có số lượng quan sát tương đối thấp. Nếu bạn muốn cải thiện lượng thông tin bạn có thể nhận được từ biến này, bạn có thể đúc lại nó ở cấp cao hơn. Cụ thể, bạn tạo các nhóm lớn hơn với trình độ học vấn tương tự. Ví dụ, trình độ học vấn thấp sẽ được chuyển đổi khi bỏ học. Các cấp học cao hơn sẽ được thay đổi thành thạc sĩ.
Đây là chi tiết:
 

 
 

Cấp độ cũ
 

Cấp độ mới
 

 

 

 
 

Trường mầm non
 

rơi ra ngoài
 

 
 

ngày 10
 

Rơi ra ngoài
 

 
 

Ngày 11
 

Rơi ra ngoài
 

 
 

Ngày 12
 

Rơi ra ngoài
 

 
 

1-4
 

Rơi ra ngoài
 

 
 

5-6
 

Rơi ra ngoài
 

 
 

7-8
 

Rơi ra ngoài
 

 
 

Thứ 9
 

Rơi ra ngoài
 

 
 

HS-Grad
 

HighGrad
 

 
 

Đại học nào đó
 

cộng đồng
 

 
 

PGS-acdm
 

cộng đồng
 

 
 

PGS-voc
 

cộng đồng
 

 
 

Cử nhân
 

Cử nhân
 

 
 

Thạc sĩ
 

Thạc sĩ
 

 
 

Trường sư phạm
 

Thạc sĩ
 

 
 

Tiến sĩ
 

Bằng tiến sĩ
 

 

 

 
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")))))))
Giải thích mã

Chúng tôi sử dụng động từ đột biến từ thư viện dplyr. Chúng tôi thay đổi các giá trị của giáo dục bằng câu lệnh ifelse

Trong bảng dưới đây, bạn tạo một thống kê tóm tắt để xem trung bình, cần bao nhiêu năm học (z-value) để đạt được bằng Cử nhân, Thạc sĩ hoặc Tiến sĩ.
recast_data % > %group_by(education) % > %summarize(average_educ_year = mean(educational.num),count = n()) % > %arrange(average_educ_year)
Đầu ra:
## # A tibble: 6 x 3## education average_educ_year count##   ## 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
Đọc lại tình trạng hôn nhân
Cũng có thể tạo các cấp độ thấp hơn cho tình trạng hôn nhân. Trong đoạn mã sau, bạn thay đổi cấp độ như sau:
 

 
 

Cấp độ cũ
 

Cấp độ mới
 

 

 

 
 

Không bao giờ kết hôn
 

Chưa kết hôn
 

 
 

Đã kết hôn-vợ / chồng-vắng mặt
 

Chưa kết hôn
 

 
 

Đã kết hôn-AF-vợ / chồng
 

Cưới nhau
 

 
 

Đã kết hôn-công dân-vợ / chồng
 
 

 
 

Ly thân
 

Ly thân
 

 
 

Đã ly hôn
 
 

 
 

Góa phụ
 

Góa phụ
 

 

 

 
# Change level marryrecast_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")))))
 Bạn có thể kiểm tra số lượng cá nhân trong mỗi nhóm. 
table(recast_data$marital.status)
Đầu ra:
## ## Married Not_married Separated Widow## 21165 15359 7727 1286 
Bước 4) Thống kê Tóm tắt
Đã đến lúc kiểm tra một số thống kê về các biến mục tiêu của chúng tôi. Trong biểu đồ bên dưới, bạn đếm tỷ lệ phần trăm cá nhân kiếm được hơn 50k dựa trên giới tính của họ.
# Plot gender incomeggplot(recast_data, aes(x = gender, fill = income)) +geom_bar(position = "fill") +theme_classic()
Đầu ra:

Tiếp theo, hãy kiểm tra xem nguồn gốc của cá nhân có ảnh hưởng đến thu nhập của họ hay không.
# Plot origin incomeggplot(recast_data, aes(x = race, fill = income)) +geom_bar(position = "fill") +theme_classic() +theme(axis.text.x = element_text(angle = 90))
Đầu ra:

Số giờ làm việc theo giới tính.
# box plot gender working timeggplot(recast_data, aes(x = gender, y = hours.per.week)) +geom_boxplot() +stat_summary(fun.y = mean,geom = "point",size = 3,color = "steelblue") +theme_classic()
Đầu ra:

Biểu đồ hộp xác nhận rằng phân phối thời gian làm việc phù hợp với các nhóm khác nhau. Trong biểu đồ hộp, cả hai giới tính không có quan sát đồng nhất.
Bạn có thể kiểm tra mật độ thời gian làm việc hàng tuần theo loại hình giáo dục. Các bản phân phối có nhiều lựa chọn riêng biệt. Nó có thể được giải thích bởi loại hợp đồng ở Mỹ.
# Plot distribution working time by educationggplot(recast_data, aes(x = hours.per.week)) +geom_density(aes(color = education), alpha = 0.5) +theme_classic()
Giải thích mã

ggplot (recast_data, aes (x = hours.per.week)): Biểu đồ mật độ chỉ yêu cầu một biến
geom_density (aes (color = education), alpha = 0,5): Đối tượng hình học để kiểm soát mật độ

Đầu ra:

Để xác nhận suy nghĩ của mình, bạn có thể thực hiện kiểm tra ANOVA một chiều:
anova <- aov(hours.per.week~education, recast_data)summary(anova)
Đầu ra:
## 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
Kiểm tra ANOVA xác nhận sự khác biệt về mức trung bình giữa các nhóm.
Không tuyến tính
Trước khi chạy mô hình, bạn có thể xem số giờ làm việc có liên quan đến tuổi tác hay không.
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()
Giải thích mã

ggplot (recast_data, aes (x = age, y = hours.per.week)): Đặt tính thẩm mỹ của biểu đồ
geom_point (aes (color = thu nhập), size = 0,5): Xây dựng biểu đồ chấm
stat_smooth (): Thêm đường xu hướng với các đối số sau:

method = 'lm': Vẽ đồ thị giá trị phù hợp nếu hồi quy tuyến tính
công thức = y ~ poly (x, 2): Phù hợp với một hồi quy đa thức
se = TRUE: Thêm lỗi chuẩn
aes (màu = thu nhập): Chia mô hình theo thu nhập



Đầu ra:

Tóm lại, bạn có thể kiểm tra các điều khoản tương tác trong mô hình để nhận ra hiệu ứng không tuyến tính giữa thời gian làm việc hàng tuần và các tính năng khác. Điều quan trọng là phát hiện thời gian làm việc khác nhau trong điều kiện nào.
Tương quan
Việc kiểm tra tiếp theo là hình dung mối tương quan giữa các biến. Bạn chuyển đổi loại mức yếu tố thành số để bạn có thể vẽ bản đồ nhiệt có chứa hệ số tương quan được tính bằng phương pháp Spearman.
library(GGally)# Convert data to numericcorr <- 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")
Giải thích mã

data.frame (lapply (recast_data, as.integer)): Chuyển đổi dữ liệu thành số
ggcorr () vẽ bản đồ nhiệt với các đối số sau:

method: Phương pháp tính toán mối tương quan
nbreaks = 6: Số lần ngắt
hjust = 0.8: Vị trí điều khiển của tên biến trong biểu đồ
label = TRUE: Thêm nhãn vào giữa cửa sổ
label_size = 3: Nhãn kích thước
color = "grey50"): Màu của nhãn



Đầu ra:

Bước 5) Đào tạo / bộ kiểm tra
Mọi tác vụ học máy được giám sát đều yêu cầu phân tách dữ liệu giữa tập hợp tàu và tập thử nghiệm. Bạn có thể sử dụng "chức năng" mà bạn đã tạo trong các hướng dẫn học tập có giám sát khác để tạo tập huấn luyện / kiểm tra.
set.seed(1234)create_train_test <- function(data, size = 0.8, train = TRUE) {n_row = nrow(data)total_row = size * n_rowtrain_sample <- 1: total_rowif (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)
Đầu ra:
## [1] 36429 9
 
dim(data_test)
Đầu ra:
## [1] 9108 9 
Bước 6) Xây dựng mô hình
Để xem thuật toán hoạt động như thế nào, bạn sử dụng gói glm (). Các tổng quát mô hình tuyến tính là một tập hợp các mô hình. Cú pháp cơ bản là:
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")
Bạn đã sẵn sàng ước tính mô hình hậu cần để phân chia mức thu nhập giữa một tập hợp các tính năng.
formula <- income~.logit <- glm(formula, data = data_train, family = 'binomial')summary(logit)
Giải thích mã

công thức <- thu nhập ~ .: Tạo mô hình để phù hợp
logit <- glm (công thức, dữ liệu = data_train, family = 'binomial'): Điều chỉnh mô hình logistic (family = 'binomial') với dữ liệu data_train.
tóm tắt (logit): In tóm tắt của mô hình

Đầu ra:
#### 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
Bản tóm tắt mô hình của chúng tôi tiết lộ thông tin thú vị. Hiệu suất của hồi quy logistic được đánh giá bằng các số liệu chính cụ thể.

AIC (Tiêu chí thông tin Akaike): Đây là tương đương với R2 trong hồi quy logistic. Nó đo lường sự phù hợp khi một hình phạt được áp dụng cho số lượng tham số. Giá trị AIC nhỏ hơn cho thấy mô hình gần với sự thật hơn.
Độ lệch rỗng: Chỉ phù hợp với mô hình với phần đánh chặn. Bậc tự do là n-1. Chúng ta có thể hiểu nó như một giá trị Chi-square (giá trị phù hợp khác với kiểm định giả thuyết giá trị thực tế).
Độ lệch dư: Mô hình với tất cả các biến. Nó cũng được hiểu là một thử nghiệm giả thuyết Chi-square.
Số lần lặp lại Fisher Scoring: Số lần lặp lại trước khi hội tụ.

Đầu ra của hàm glm () được lưu trữ trong một danh sách. Đoạn mã dưới đây hiển thị tất cả các mục có sẵn trong biến logit mà chúng tôi đã xây dựng để đánh giá hồi quy logistic.
# Danh sách rất dài, chỉ in ba phần tử đầu tiên
lapply(logit, class)[1:3]
Đầu ra:
## $coefficients## [1] "numeric"#### $residuals## [1] "numeric"#### $fitted.values## [1] "numeric"
Mỗi giá trị có thể được trích xuất với dấu $ theo sau tên của số liệu. Ví dụ, bạn đã lưu trữ mô hình dưới dạng logit. Để trích xuất các tiêu chí AIC, bạn sử dụng:
logit$aic
Đầu ra:
## [1] 27086.65
Bước 7) Đánh giá hiệu suất của mô hình
Ma trận hỗn loạn
Các ma trận nhầm lẫn là một lựa chọn tốt hơn để đánh giá việc thực hiện phân loại so với các số liệu khác nhau mà bạn thấy trước đó. Ý tưởng chung là đếm số lần các trường hợp Đúng được phân loại là Sai.

Để tính toán ma trận nhầm lẫn, trước tiên bạn cần có một tập hợp các dự đoán để chúng có thể được so sánh với các mục tiêu thực tế.
predict <- predict(logit, data_test, type = 'response')# confusion matrixtable_mat <- table(data_test$income, predict > 0.5)table_mat
Giải thích mã

dự đoán (logit, data_test, type = 'response'): Tính toán dự đoán trên tập kiểm tra. Đặt type = 'response' để tính xác suất phản hồi.
table (data_test $ thu nhập, dự đoán> 0.5): Tính toán ma trận nhầm lẫn. dự đoán> 0,5 có nghĩa là nó trả về 1 nếu xác suất dự đoán trên 0,5, còn lại là 0.

Đầu ra:
#### FALSE TRUE## <=50K 6310 495## >50K 1074 1229
Mỗi hàng trong ma trận nhầm lẫn đại diện cho một mục tiêu thực tế, trong khi mỗi cột đại diện cho một mục tiêu dự đoán. Hàng đầu tiên của ma trận này xem xét thu nhập thấp hơn 50k (hạng Sai): 6241 được phân loại chính xác là cá nhân có thu nhập thấp hơn 50k ( Đúng âm ), trong khi một người còn lại được phân loại sai là trên 50k ( Sai dương tính ). Hàng thứ hai xem xét thu nhập trên 50k, loại tích cực là 1229 ( Đúng dương ), trong khi âm Đúng là 1074.
Bạn có thể tính toán độ chính xác của mô hình bằng cách tính tổng số dương thực + âm thực trên tổng số quan sát

accuracy_Test <- sum(diag(table_mat)) / sum(table_mat)accuracy_Test
Giải thích mã

sum (Diag (table_mat)): Tổng của đường chéo
sum (table_mat): Tổng của ma trận.

Đầu ra:
## [1] 0.8277339 
Mô hình dường như mắc phải một vấn đề, nó đánh giá quá cao số lượng âm tính sai. Đây được gọi là nghịch lý kiểm tra độ chính xác . Chúng tôi đã nói rằng độ chính xác là tỷ lệ các dự đoán đúng trên tổng số các trường hợp. Chúng ta có thể có độ chính xác tương đối cao nhưng là một mô hình vô dụng. Nó xảy ra khi có một giai cấp thống trị. Nếu bạn nhìn lại ma trận nhầm lẫn, bạn có thể thấy hầu hết các trường hợp được phân loại là tiêu cực thực sự. Hãy tưởng tượng bây giờ, mô hình phân loại tất cả các lớp là âm (tức là thấp hơn 50k). Bạn sẽ có độ chính xác là 75 phần trăm (6718/6718 + 2257). Mô hình của bạn hoạt động tốt hơn nhưng lại gặp khó khăn trong việc phân biệt điều tích cực thực sự với điều tiêu cực thực sự.
Trong tình huống như vậy, tốt hơn là có một số liệu ngắn gọn hơn. Chúng ta có thể xem xét:

Độ chính xác = TP / (TP + FP)
Nhớ lại = TP / (TP + FN)

Độ chính xác so với Nhớ lại
Precision xem xét độ chính xác của dự đoán tích cực. Nhớ lại là tỷ lệ các trường hợp tích cực được phát hiện chính xác bởi bộ phân loại;
Bạn có thể xây dựng hai hàm để tính toán hai số liệu này

Xây dựng độ chính xác

precision <- function(matrix) {# True positivetp <- matrix[2, 2]# false positivefp <- matrix[1, 2]return (tp / (tp + fp))}
Giải thích mã

mat [1,1]: Trả về ô đầu tiên của cột đầu tiên của khung dữ liệu, tức là giá trị dương thực sự
mat [1,2]; Trả về ô đầu tiên của cột thứ hai của khung dữ liệu, tức là giá trị dương giả

recall <- function(matrix) {# true positivetp <- matrix[2, 2]# false positivefn <- matrix[2, 1]return (tp / (tp + fn))}
Giải thích mã

mat [1,1]: Trả về ô đầu tiên của cột đầu tiên của khung dữ liệu, tức là giá trị dương thực sự
mat [2,1]; Trả về ô thứ hai của cột đầu tiên của khung dữ liệu, tức là giá trị âm sai

Bạn có thể kiểm tra các chức năng của mình
prec <- precision(table_mat)precrec <- recall(table_mat)rec
Đầu ra:
## [1] 0.712877## [2] 0.5336518
Khi mô hình cho biết đó là một cá nhân trên 50k, nó chỉ đúng trong 54% trường hợp và có thể xác nhận những cá nhân trên 50k trong 72% trường hợp.
Bạn có thể tạo trị trung bình hài hòa của hai chỉ số này, có nghĩa là nó có trọng số nhiều hơn cho các giá trị thấp hơn.

f1 <- 2 * ((prec * rec) / (prec + rec))f1
Đầu ra:
## [1] 0.6103799 
Sự cân bằng giữa chính xác và thu hồi
Không thể có cả độ chính xác cao và thu hồi cao.
Nếu chúng tôi tăng độ chính xác, cá nhân chính xác sẽ được dự đoán tốt hơn, nhưng chúng tôi sẽ bỏ lỡ rất nhiều trong số đó (thu hồi thấp hơn). Trong một số tình huống, chúng tôi thích độ chính xác cao hơn là thu hồi. Có một mối quan hệ lõm giữa độ chính xác và thu hồi.

Hãy tưởng tượng, bạn cần dự đoán xem một bệnh nhân có mắc bệnh hay không. Bạn muốn chính xác nhất có thể.
Nếu bạn cần phát hiện những kẻ gian lận tiềm năng trên đường phố thông qua nhận dạng khuôn mặt, tốt hơn hết là bạn nên bắt nhiều người bị gắn mác lừa đảo mặc dù độ chính xác thấp. Cảnh sát sẽ có thể thả cá nhân không gian lận.

Đường cong ROC
Các hoạt động tiếp nhận đường cong là một công cụ thường được sử dụng với phân loại nhị phân. Nó rất giống với đường cong độ chính xác / thu hồi, nhưng thay vì vẽ biểu đồ độ chính xác so với thu hồi, đường cong ROC hiển thị tỷ lệ dương tính thực sự (tức là thu hồi) so với tỷ lệ dương tính giả. Tỷ lệ dương tính giả là tỷ lệ các trường hợp tiêu cực được phân loại không chính xác là dương tính. Nó bằng một trừ đi tỷ lệ âm thực sự. Tỷ lệ âm thực sự còn được gọi là độ đặc hiệu . Do đó, đường cong ROC biểu thị độ nhạy (thu hồi) so với độ đặc hiệu 1
Để vẽ đường cong ROC, chúng ta cần cài đặt một thư viện có tên là RORC. Chúng tôi có thể tìm thấy trong thư viện chung cư. Bạn có thể nhập mã:
conda install -cr r-rocr --yes
Chúng ta có thể vẽ biểu đồ ROC bằng các hàm dự đoán () và hiệu suất ().
library(ROCR)ROCRpred <- prediction(predict, data_test$income)ROCRperf <- performance(ROCRpred, 'tpr', 'fpr')plot(ROCRperf, colorize = TRUE, text.adj = c(-0.2, 1.7))
Giải thích mã

dự đoán (dự đoán, data_test $ thu nhập): Thư viện ROCR cần tạo một đối tượng dự đoán để chuyển đổi dữ liệu đầu vào
hiệu suất (ROCRpred, 'tpr', 'fpr'): Trả về hai kết hợp để tạo ra trong biểu đồ. Ở đây, tpr và fpr được xây dựng. Toàn bộ độ chính xác của đồ thị và gọi lại cùng nhau, sử dụng "Pre", "rec".

Đầu ra:

Bước 8) Cải thiện mô hình
Bạn có thể cố gắng thêm tính không tuyến tính vào mô hình với sự tương tác giữa

tuổi và giờ.per.week
giới tính và giờ.per.week.

Bạn cần sử dụng bài kiểm tra điểm số để so sánh cả hai mô hình
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
Đầu ra:
## [1] 0.6109181 
Điểm số cao hơn một chút so với lần trước. Bạn có thể tiếp tục làm việc trên dữ liệu để cố gắng vượt qua điểm số.

Tóm lược
Chúng ta có thể tóm tắt hàm để huấn luyện hồi quy logistic trong bảng dưới đây:
 

 
 

Gói hàng
 

Mục tiêu
 

chức năng
 

tranh luận
 

 

 

 
 

-
 

Tạo tập dữ liệu huấn luyện / kiểm tra
 

create_train_set ()
 

dữ liệu, kích thước, chuyến tàu
 

 
 

glm
 

Đào tạo một mô hình tuyến tính tổng quát
 

glm ()
 

công thức, dữ liệu, họ *
 

 
 

glm
 

Tóm tắt mô hình
 

tóm lược()
 

mô hình vừa vặn
 

 
 

căn cứ
 

Đưa ra dự đoán
 

dự đoán ()
 

mô hình được trang bị, tập dữ liệu, type = 'response'
 

 
 

căn cứ
 

Tạo một ma trận nhầm lẫn
 

bàn()
 

y, dự đoán ()
 

 
 

căn cứ
 

Tạo điểm chính xác
 

sum (đường chéo (table ()) / sum (table ()
 

 
 

ROCR
 

Tạo ROC: Bước 1 Tạo dự đoán
 

sự dự đoán()
 

dự đoán (), y
 

 
 

ROCR
 

Tạo ROC: Bước 2 Tạo hiệu suất
 

hiệu suất()
 

dự đoán (), 'tpr', 'fpr'
 

 
 

ROCR
 

Tạo ROC: Bước 3 Vẽ biểu đồ
 

âm mưu()
 

hiệu suất()
 

 

 

Các loại mô hình GLM khác là:
- nhị thức: (link = "logit")
- gaussian: (link = "danh tính")
- Gamma: (link = "inverse")
- inverse.gaussian: (link = "1 / mu 2")
- poisson: (link = "log")
- gần như: (liên kết = "nhận dạng", phương sai = "hằng số")
- quasibinomial: (link = "logit")
- quasipoisson: (link = "log")