Nhiệm vụ tuần 1

Nhiệm vụ tuần 1

Tóm tắt Sách “Generalized Linear Models With Examples in R”

Giới thiệu Chung

Cuốn sách này trình bày về Các Mô hình Tuyến tính Tổng quát (Generalized Linear Models - GLM) với các ví dụ minh họa sử dụng ngôn ngữ lập trình R [1]. Sách được chuẩn bị bằng LATEX và R phiên bản 3.4.3, tích hợp sử dụng Sweave [2]. Cuốn sách này cung cấp dữ liệu qua gói R GLMsData [1, 2]. Các hàm trong R được hiển thị bằng phông chữ typewriter theo sau bởi dấu ngoặc đơn, ví dụ: glm() [1]. Các toán tử, data frame và biến trong R cũng được hiển thị bằng phông chữ typewriter, ví dụ: Smoke [1]. Các gói R được hiển thị bằng phông chữ in đậm sans serif, ví dụ: GLMsData [1].

Chương 1: Mô hình Thống kê (Statistical Models)

Chương này giới thiệu các khái niệm cơ bản về mô hình thống kê [3].

• Quy ước mô tả dữ liệu: Giới thiệu cách mô tả dữ liệu, thường sử dụng ví dụ minh họa [4]. Ví dụ 1.1 sử dụng bộ dữ liệu lungcap từ gói GLMsData, nghiên cứu mối quan hệ giữa dung tích phổi (fev), tình trạng hút thuốc, tuổi, chiều cao và giới tính trên 654 thanh niên ở Đông Boston [4, 5]. Có thể nạp dữ liệu này vào R bằng cách tải gói GLMsData và sử dụng lệnh data(lungcap) [4]. Có thể xem các giá trị đầu tiên/cuối cùng của một biến bằng head() và tail() [5].

• Vẽ đồ thị dữ liệu: Trình bày cách vẽ đồ thị để khám phá dữ liệu [3, 6-11]. R có các cơ chế mạnh mẽ để tạo đồ thị [12]. Các đồ thị đơn giản có thể được tạo dễ dàng [12]. Có thể tùy chỉnh đồ thị bằng các tham số đồ họa của R như las, ylim, xlim, xlab, ylab, main, pch [13]. Có thể sử dụng lệnh plot(lungcap$FEV ~ lungcap$Age) để vẽ đồ thị FEV so với Tuổi [12].

• Mô hình thống kê mô tả cả đặc điểm ngẫu nhiên và có hệ thống của dữ liệu [14].

• Mô hình hồi quy (Regression Models) [14].

• Giải thích mô hình hồi quy [14].

• Tất cả mô hình đều sai, nhưng một số hữu ích [14].

• Mục đích của mô hình thống kê ảnh hưởng đến cách phát triển nó [14].

• Độ chính xác so với Tính tiết kiệm (Accuracy vs Parsimony) [14].

• Thí nghiệm so với Nghiên cứu quan sát: Nguyên nhân so với Liên kết (Causality vs Association) [15].

• Thu thập dữ liệu và Khả năng tổng quát hóa [15].

• Sử dụng R cho mô hình thống kê [15].

Chương 2: Mô hình Hồi quy Tuyến tính (Linear Regression Models)

Chương này đi sâu vào mô hình hồi quy tuyến tính [16].

• Định nghĩa mô hình hồi quy tuyến tính [16].

• Hồi quy tuyến tính đơn giản [16]. Bao gồm Ước lượng Bình phương Tối thiểu (Least-Squares Estimation), Ước lượng Hệ số, Ước lượng Phương sai sigma^2, Sai số Chuẩn của Hệ số, Sai số Chuẩn của Giá trị Fitted [16].

• Ước lượng cho hồi quy bội (Multiple Regression) [17]. Bao gồm Ước lượng Hệ số, Ước lượng Phương sai sigma^2, Sai số Chuẩn [17].

• Công thức Ma trận của Mô hình Hồi quy Tuyến tính [17]. Bao gồm Ký hiệu Ma trận, Ước lượng Hệ số, Ước lượng Phương sai sigma^2, Ước lượng Phương sai của beta-hat, Ước lượng Phương sai của Giá trị Fitted [17]. Ví dụ 2.10 cho thấy cách sử dụng R để tạo ma trận thiết kế Xmat bằng model.matrix(), tính toán X^T X và X^T y sử dụng t() (chuyển vị) và %*% (nhân ma trận), và tìm ước lượng hệ số beta bằng cách giải hệ phương trình sử dụng solve() (nghịch đảo ma trận) [18].

• Phân tích Phương sai cho Mô hình Hồi quy [19].

• So sánh các Mô hình Lồng nhau (Nested Models) [19]. Có thể sử dụng phân tích phương sai để so sánh hai mô hình lồng nhau [19]. R lệnh anova() có thể được sử dụng cho mục đích này [20, 21]. Bảng ANOVA bao gồm các cột Df (Bậc tự do), Sum Sq (Tổng bình phương), Mean Sq (Bình phương trung bình), F value (Giá trị F) và Pr(>F) (P-value) [21].

• Chẩn đoán cho Mô hình Hồi quy Tuyến tính [22].

• Xây dựng Mô hình Tuyến tính [22].

• Sử dụng R cho Phân tích Chẩn đoán của Mô hình Hồi quy Tuyến tính [22].

• Hệ số xác định R-squared (R^2): Tỷ lệ biến thiên tổng cộng được giải thích bởi hồi quy [20]. Có thể tính R^2 bằng công thức 1 - RSS/SST hoặc lấy trực tiếp từ summary(model)$r.squared [23]. R^2 hiệu chỉnh (adjusted R^2) cũng được báo cáo trong summary() [23].

• So sánh các mô hình không lồng nhau có thể sử dụng Tiêu chí Thông tin Akaike (AIC) hoặc Tiêu chí Thông tin Bayes (BIC) [24]. AIC có thể được trích xuất bằng lệnh extractAIC() [25]. Đối với mô hình hồi quy tuyến tính, bậc tự do tương đương là số lượng tham số hồi quy được ước lượng [25]. BIC sử dụng cùng hàm extractAIC() nhưng với hình phạt (k) được điều chỉnh [25].

• Problem 2.14 mô tả một nghiên cứu về huyết áp ở nam giới Ghana, sử dụng mô hình hồi quy tuyến tính với nhiều biến giải thích và yêu cầu so sánh các mô hình bằng ANOVA, AIC, BIC, R^2 và R^2 hiệu chỉnh [24].

Chương 3 & 4: Chẩn đoán và Xây dựng Mô hình Hồi quy Tuyến tính

Các chương này mở rộng về chẩn đoán và xây dựng mô hình [22, 26].

• Chẩn đoán mô hình bao gồm kiểm tra các giả định [22, 27-29]. Các đồ thị chẩn đoán như phần dư chuẩn hóa (standardized residuals) so với giá trị fitted là hữu ích [27-29].

• Biến đổi đơn giản các biến giải thích (Simple Transformations of Covariates) [30]. Có thể biến đổi biến giải thích, ví dụ sử dụng logarit (log(Ht)) [30].

• Biến đổi biến phản hồi (Transformations of the Response) [27-29]. Có thể kiểm tra tác động của các biến đổi đối với biến phản hồi bằng cách xem xét các đồ thị phần dư [27-29].

• Nguyên tắc biên (Marginality Principle) [19, 31].

• Phân tích ảnh hưởng (Influence Analysis) [32]. Có thể sử dụng influence.measures() và cooks.distance() trong R để xác định các quan sát có ảnh hưởng lớn [32, 33]. Quan sát có đòn bẩy (leverage) cao có thể không có ảnh hưởng nếu phần dư nhỏ, và ngược lại [33]. Có thể cập nhật mô hình bằng cách loại bỏ các quan sát có ảnh hưởng lớn bằng tùy chọn subset trong hàm mô hình [32].

• Lựa chọn biến (Variable Selection) [19].

• Problem 3.8 sử dụng dữ liệu lungcap để so sánh các mô hình hồi quy tuyến tính khác nhau [34].

• Problem 3.23 sử dụng dữ liệu gopher để khám phá mối quan hệ giữa các biến, fit mô hình hồi quy tuyến tính có trọng số (weights), và thực hiện phân tích chẩn đoán [35].

• Problem 3.25 sử dụng dữ liệu ratliver để khám phá dữ liệu, fit mô hình hồi quy tuyến tính, và kiểm tra các quan sát có ảnh hưởng [32].

Chương 5: Mô hình Tuyến tính Tổng quát: Cấu trúc (Generalized Linear Models: Structure)

Chương này giới thiệu cấu trúc của GLM [26].

• Hai thành phần của Mô hình Tuyến tính Tổng quát: Thành phần ngẫu nhiên (Random Component) và thành phần có hệ thống (Systematic Component) liên kết thông qua hàm liên kết (link function) [26, 36].

• Thành phần ngẫu nhiên: Mô hình Phân tán Hàm mũ (Exponential Dispersion Models - EDM) [26, 36]. Bảng 5.1 liệt kê các EDM phổ biến bao gồm Normal, Binomial, Negative Binomial, Poisson, Gamma, Inverse Gaussian và Tweedie, cùng với các thuộc tính của chúng như hàm phương sai V(μ), hàm tích lũy κ(θ), tham số chính tắc θ, tham số phân tán φ, độ lệch đơn vị d(y,μ), miền hỗ trợ S và miền cho μ và θ [36-38]. Phân phối Tweedie bao gồm các trường hợp đặc biệt như Gamma (ξ=2) và Poisson (ξ=1 với φ=1) [36].

Các Chương Khác về GLM Chuyên biệt

Cuốn sách tiếp tục khám phá các loại GLM cụ thể.

• Chương 9: Models for Proportions: Binomial GLMs [37]. Sử dụng cho dữ liệu tỷ lệ [39]. Ví dụ sử dụng dữ liệu germ với mô hình Binomial GLM có trọng số (weights) [39].

• Chương 10: Models for Counts: Poisson and Negative Binomial GLMs [38, 40]. Sử dụng cho dữ liệu đếm [40, 41].

  • Poisson GLM: Thích hợp cho dữ liệu đếm khi phương sai bằng trung bình [40]. Có thể bao gồm một offset trong mô hình, ví dụ:offset(log(Pop)) trong mô hình Poisson log μ = log T + ... nơi log T là một offset không cần ước lượng tham số [40]. Ví dụ sử dụng dữ liệu danishlc [40].

  • Negative Binomial GLM: Thích hợp cho dữ liệu đếm khi có hiện tượng overdispersion (phương sai lớn hơn trung bình) [42]. Có thể so sánh mô hình Quasi-poisson và Negative Binomial khi có overdispersion [42].

  • Quasi-poisson GLM: Cũng được sử dụng khi có overdispersion trong dữ liệu đếm [43]. anova() với test="F" thường được sử dụng cho mô hình Quasi-poisson [43]. Overdispersion có thể được xác nhận nếu độ lệch dư (residual deviance) hoặc Pearson Chi-squared lớn hơn đáng kể bậc tự do dư (residual degrees of freedom) [43]. Ví dụ sử dụng dữ liệu hcrabs [43].

  • Kiểm định độ khớp (Goodness-of-fit): Có thể sử dụng pchisq(deviance(model), df.residual(model), lower.tail=FALSE) để kiểm định giả thuyết rằng mô hình là phù hợp [41].

• Chương 11: Positive Continuous Data: Gamma and Inverse Gaussian GLMs [38, 44]. Sử dụng cho dữ liệu liên tục dương [44]. Ví dụ sử dụng dữ liệu lime cho Gamma GLM [44, 45].

• Chương 12: Tweedie GLMs [46]. Trình bày các Mô hình Phân tán Hàm mũ Tweedie (Tweedie EDMs) và cách fit Mô hình Tuyến tính Tổng quát Tweedie (Tweedie GLMs) [38, 46, 47]. Bao gồm ước lượng tham số chỉ số ξ [47]. Có thể sử dụng hàm tweedie() và tweedie.profile() trong R [48]. summary() cho mô hình Tweedie sẽ hiển thị ước lượng, sai số chuẩn, giá trị t và P-value, cùng với ước lượng tham số phân tán [49].

Suy luận cho GLM (Inference for GLMs)

• Kiểm định Tỷ số Độ Lệch (Deviance Ratio Test) hay Kiểm định Tỷ số Khả Năng (Likelihood Ratio Test) [50, 51]. Đối với các phân phối có φ = 1 (ví dụ: Poisson), sự khác biệt giữa độ lệch dư của hai mô hình lồng nhau theo phân phối chi-squared [50, 51]. Có thể tính toán sự khác biệt này và P-value bằng deviance() và pchisq() [50]. Lệnh anova(model, test="Chisq") trong R thực hiện kiểm định này [51].

• Kiểm định Wald (Wald Test): Sử dụng cho các tham số mô hình riêng lẻ hoặc các giả thuyết cụ thể [52, 53]. R summary() cho GLM thường thực hiện z-test (coi phân tán đã biết), điều này có thể phóng đại nhẹ ý nghĩa thống kê [52]. Có thể tính toán kiểm định Wald thủ công [52].

• Kiểm định Score (Score Test): Một phương pháp kiểm định khác, có thể thực hiện trong R bằng hàm như glm.scoretest từ gói statmod [54].

Phụ lục A: Sử dụng R cho Phân tích Dữ liệu (Using R for Data Analysis)

Phụ lục này giới thiệu cách sử dụng R [55].

• Chuẩn bị sử dụng R: Cách tải và cài đặt R và các gói R cần thiết (ví dụ: GLMsData) [4, 55, 56]. Các trang web R quan trọng bao gồm CRAN và R-project.org [56-58].

• Giới thiệu cơ bản về sử dụng R: Sử dụng R như một máy tính nâng cao, thoát R (q()), nhận trợ giúp (help()), tên biến, làm việc với vector, nạp dữ liệu vào R [55, 59]. Có thể nạp dữ liệu có sẵn trong R bằng lệnh data() [60]. Có thể nạp dữ liệu từ file bằng read.table() [61].

• Làm việc với data frame: Ví dụ lungcap là một data frame [4]. Có thể truy cập biến trong data frame bằng $, ví dụ lungcap$Age [12, 61]. Có thể sử dụng with() để truy cập biến dễ dàng hơn, nhưng attach()/detach() không được khuyến khích [61]. Có thể xem cấu trúc data frame bằng str(), tóm tắt bằng summary() [48].

• Sử dụng hàm trong R: R chứa rất nhiều hàm và các gói bổ sung thêm nhiều hàm khác [59, 61]. Nhiều hàm đã được sử dụng trong sách như glm(), lm(), anova(), summary(), plot() [1, 20, 48].

• Các hàm thống kê cơ bản trong R: mean(), sd(), var(), median(), IQR() [12, 59, 61].

• Vẽ đồ thị cơ bản trong R [12, 13, 59, 62-64].

• Viết hàm trong R: Hướng dẫn cách viết hàm tùy chỉnh, sử dụng đối số mặc định, lệnh signif(), paste(), return() [59, 65, 66].

Phụ lục B: Bộ dữ liệu (Data sets)

Liệt kê các bộ dữ liệu được sử dụng trong sách, bao gồm lungcap, hcrabs, lime, pock, gopher, ratliver, blocks, danishlc, cervical, turbines,… [67].

Chỉ mục (Index)

Bao gồm chỉ mục các bộ dữ liệu và các lệnh/hàm trong R được sử dụng trong sách, cung cấp tham chiếu trang [48, 68, 69].

Các Bài toán bổ sung (Extra Problems)

Cung cấp các bài toán thực hành bổ sung [70, 71].

Thực hiện thống kê mô tả cho các biến trong file: Supermarket Transactions.csv

if (!require("readxl")) install.packages("readxl")

## Loading required package: readxl

if (!require("dplyr")) install.packages("dplyr")

## Loading required package: dplyr

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

if (!require("psych")) install.packages("psych")

## Loading required package: psych

  library(readxl)
data <- read_excel("C:/Users/Nep/Downloads/Supermarket Transactions.xlsx")

## New names:
## • `` -> `...1`

str(data)

## tibble [14,059 × 16] (S3: tbl_df/tbl/data.frame)
##  $ ...1             : num [1:14059] 1 2 3 4 5 6 7 8 9 10 ...
##  $ PurchaseDate     : POSIXct[1:14059], format: "2007-12-18" "2007-12-20" ...
##  $ CustomerID       : num [1:14059] 7223 7841 8374 9619 1900 ...
##  $ Gender           : chr [1:14059] "F" "M" "F" "M" ...
##  $ MaritalStatus    : chr [1:14059] "S" "M" "M" "M" ...
##  $ Homeowner        : chr [1:14059] "Y" "Y" "N" "Y" ...
##  $ Children         : num [1:14059] 2 5 2 3 3 3 2 2 3 1 ...
##  $ AnnualIncome     : chr [1:14059] "$30K - $50K" "$70K - $90K" "$50K - $70K" "$30K - $50K" ...
##  $ City             : chr [1:14059] "Los Angeles" "Los Angeles" "Bremerton" "Portland" ...
##  $ StateorProvince  : chr [1:14059] "CA" "CA" "WA" "OR" ...
##  $ Country          : chr [1:14059] "USA" "USA" "USA" "USA" ...
##  $ ProductFamily    : chr [1:14059] "Food" "Food" "Food" "Food" ...
##  $ ProductDepartment: chr [1:14059] "Snack Foods" "Produce" "Snack Foods" "Snacks" ...
##  $ ProductCategory  : chr [1:14059] "Snack Foods" "Vegetables" "Snack Foods" "Candy" ...
##  $ UnitsSold        : num [1:14059] 5 5 3 4 4 3 4 6 1 2 ...
##  $ Revenue          : num [1:14059] 27.38 14.9 5.52 4.44 14 ...

head(data)

numeric_vars <- data[, sapply(data, is.numeric)]
summary(numeric_vars)

##       ...1         CustomerID       Children      UnitsSold        Revenue     
##  Min.   :    1   Min.   :    3   Min.   :0.00   Min.   :1.000   Min.   : 0.53  
##  1st Qu.: 3516   1st Qu.: 2549   1st Qu.:1.00   1st Qu.:3.000   1st Qu.: 6.84  
##  Median : 7030   Median : 5060   Median :3.00   Median :4.000   Median :11.25  
##  Mean   : 7030   Mean   : 5117   Mean   :2.53   Mean   :4.081   Mean   :13.00  
##  3rd Qu.:10544   3rd Qu.: 7633   3rd Qu.:4.00   3rd Qu.:5.000   3rd Qu.:17.37  
##  Max.   :14059   Max.   :10280   Max.   :5.00   Max.   :8.000   Max.   :56.70

library(psych)
describe(numeric_vars)

categorical_vars <- data[, sapply(data, is.character) | sapply(data, is.factor)]
lapply(categorical_vars, table)

## $Gender
## 
##    F    M 
## 7170 6889 
## 
## $MaritalStatus
## 
##    M    S 
## 6866 7193 
## 
## $Homeowner
## 
##    N    Y 
## 5615 8444 
## 
## $AnnualIncome
## 
##   $10K - $30K $110K - $130K $130K - $150K       $150K +   $30K - $50K 
##          3090           643           760           273          4601 
##   $50K - $70K   $70K - $90K  $90K - $110K 
##          2370          1709           613 
## 
## $City
## 
##      Acapulco    Bellingham Beverly Hills     Bremerton       Camacho 
##           383           143           811           834           452 
##   Guadalajara       Hidalgo   Los Angeles        Merida   Mexico City 
##            75           845           926           654           194 
##       Orizaba      Portland         Salem    San Andres     San Diego 
##           464           876          1386           621           866 
## San Francisco       Seattle       Spokane        Tacoma     Vancouver 
##           130           922           875          1257           633 
##      Victoria   Walla Walla        Yakima 
##           176           160           376 
## 
## $StateorProvince
## 
##        BC        CA        DF  Guerrero   Jalisco        OR  Veracruz        WA 
##       809      2733       815       383        75      2262       464      4567 
##   Yucatan Zacatecas 
##       654      1297 
## 
## $Country
## 
## Canada Mexico    USA 
##    809   3688   9562 
## 
## $ProductFamily
## 
##          Drink           Food Non-Consumable 
##           1250          10153           2656 
## 
## $ProductDepartment
## 
## Alcoholic Beverages         Baked Goods        Baking Goods           Beverages 
##                 356                 425                1072                 680 
##     Breakfast Foods        Canned Foods     Canned Products            Carousel 
##                 188                 977                 109                  59 
##            Checkout               Dairy                Deli                Eggs 
##                  82                 903                 699                 198 
##        Frozen Foods  Health and Hygiene           Household                Meat 
##                1382                 893                1420                  89 
##         Periodicals             Produce             Seafood         Snack Foods 
##                 202                1994                 102                1600 
##              Snacks       Starchy Foods 
##                 352                 277 
## 
## $ProductCategory
## 
##         Baking Goods    Bathroom Products        Beer and Wine 
##                  484                  365                  356 
##                Bread      Breakfast Foods              Candles 
##                  425                  417                   45 
##                Candy     Canned Anchovies         Canned Clams 
##                  352                   44                   53 
##       Canned Oysters      Canned Sardines        Canned Shrimp 
##                   35                   40                   38 
##          Canned Soup          Canned Tuna Carbonated Beverages 
##                  404                   87                  154 
##    Cleaning Supplies        Cold Remedies                Dairy 
##                  189                   93                  903 
##        Decongestants               Drinks                 Eggs 
##                   85                  135                  198 
##           Electrical      Frozen Desserts       Frozen Entrees 
##                  355                  323                  118 
##                Fruit             Hardware        Hot Beverages 
##                  765                  129                  226 
##              Hygiene     Jams and Jellies     Kitchen Products 
##                  197                  588                  217 
##            Magazines                 Meat        Miscellaneous 
##                  202                  761                   42 
##  Packaged Vegetables       Pain Relievers       Paper Products 
##                   48                  192                  345 
##                Pizza     Plastic Products Pure Juice Beverages 
##                  194                  141                  165 
##              Seafood          Side Dishes          Snack Foods 
##                  102                  153                 1600 
##            Specialty        Starchy Foods           Vegetables 
##                  289                  277                 1728