BÁO CÁO TIỂU LUẬN MÔN NGÔN NGỮ LẬP TRÌNH TRONG PHÂN TÍCH DỮ LIỆU

Giới thiệu bộ dữ liệu

Bộ dữ liệu US Accidents (2022) là tập dữ liệu quy mô lớn về tai nạn giao thông tại Hoa Kỳ, do Sobhan Moosavi (ĐH Ohio State) công bố trên Kaggle. Dữ liệu được thu thập tự động từ nhiều nguồn như Bing Maps Traffic API, MapQuest, TomTom, INRIX, Here Traffic cùng hệ thống cảm biến, camera và dữ liệu thời tiết từ NOAA và USDOT.

Mỗi bản ghi mô tả chi tiết một vụ tai nạn (thời gian, vị trí, điều kiện thời tiết, mức độ nghiêm trọng, v.v.) tại 49 bang của Hoa Kỳ.

Với tính toàn diện và đa chiều, bộ dữ liệu này là nguồn thông tin quan trọng cho phân tích thống kê, mô hình hóa rủi ro, và nghiên cứu các yếu tố ảnh hưởng đến mức độ nghiêm trọng của tai nạn giao thông.

Giới thiệu và cấu trúc bộ dữ liệu

Đọc bộ dữ liệu

Accidents <- read_csv("~/US_Accidents_2022.csv")

Câu lệnh: (1) câu lệnh thuộc gói readr được dùng để đọc tệp dữ liệu định dạng CSV.

Xem cấu trúc bộ dữ liệu

glimpse(Accidents )

## Rows: 1.762.452
## Columns: 17
## $ ID                <chr> "A-512230", "A-512231", "A-512232", "A-512…
## $ Source            <chr> "Source2", "Source2", "Source2", "Source2"…
## $ Severity          <dbl> 1, 1, 1, 1, 2, 1, 2, 2, 1, 1, 1, 3, 1, 1, …
## $ Start_Time        <dttm> 2022-09-08 05:49:30, 2022-09-08 02:02:05,…
## $ End_Time          <dttm> 2022-09-08 06:34:53, 2022-09-08 04:31:32,…
## $ Start_Lat         <dbl> 41,94680, 34,52117, 37,54284, 40,89663, 41…
## $ Start_Lng         <dbl> -88,20809, -117,95808, -77,44178, -81,1784…
## $ `Distance(mi)`    <dbl> 0,00, 0,00, 0,00, 0,00, 1,91, 1,59, 0,00, …
## $ Street            <chr> "Army Trail Rd", "Pearblossom Hwy", "N 2nd…
## $ City              <chr> "Bartlett", "Littlerock", "Richmond", "All…
## $ State             <chr> "IL", "CA", "VA", "OH", "OH", "PA", "OH", …
## $ Weather_Timestamp <dttm> 2022-09-08 05:52:00, 2022-09-08 01:53:00,…
## $ `Temperature(F)`  <dbl> 58, 86, 68, 62, 63, 65, 53, 58, 70, 70, 70…
## $ `Humidity(%)`     <dbl> 90, 28, 96, 86, 87, 93, 99, 97, 100, 97, 1…
## $ `Pressure(in)`    <dbl> 29,24, 27,35, 29,71, 28,71, 29,37, 29,45, …
## $ `Visibility(mi)`  <dbl> 10, 10, 10, 7, 7, 9, 7, 10, 10, 10, 10, 10…
## $ Weather_Condition <chr> "Fair", "Fair", "Mostly Cloudy", "Mostly C…

Câu lệnh: (1) được thực hiện kiểm toán nhanh cấu trúc: in số dòng/cột, kiểu dữ liệu của từng biến và một số giá trị mẫu.

Kết quả: Kết quả cho thấy bộ dữ liệu bao gồm 1.762.452 bản ghi về các vụ tai nạn giao thông tại Hoa Kỳ, bản ghi này lưu trữ 17 biến phản ánh đầy đủ ba trục thông tin cốt lõi: thời gian, không gian và thời tiết.

Quan sát nhanh bộ dữ liệu

Xem nhanh 5 dòng đầu bộ dữ liệu

head(Accidents, 5)

Câu lệnh: (1) trong R được dùng để hiển thị 5 dòng đầu tiên của bộ dữ liệu để bạn kiểm tra cấu trúc, định dạng và nội dung của các biến.

Kết quả: Cho thấy dữ liệu được đọc đúng: đủ 17 cột, mỗi dòng tương ứng một vụ tai nạn với ID, vị trí, thời gian, Severity và bối cảnh thời tiết. Đây là ảnh “soi nhanh” để xác nhận cấu trúc trước khi phân tích sâu hơn.

Xem nhanh 5 dòng cuối bộ dữ liệu

tail(Accidents, 5)

Câu lệnh: (1) trong R được dùng để hiển thị 5 dòng cuối cùng của bộ dữ liệu, giúp xem các bản ghi mới nhất hoặc kết thúc của tập dữ liệu.

Kết quả: Kết quả hiển thị 5 dòng cuối của bộ dữ liệu, giúp quan sát nhanh cấu trúc và định dạng của các biến. Các biến như ID, Source, Severity, Start_Time,… đều được lưu đúng kiểu dữ liệu.

Kích thước bảng dữ liệu

nrow(Accidents); ncol(Accidents); dim(Accidents)

## [1] 1762452

## [1] 17

## [1] 1762452      17

Câu lệnh: (1) lần lượt trả về số dòng, số cột, và kích thước đầy đủ (dòng × cột) của bảng Accidents.

Kết quả: Kết quả cho thấy tập dữ liệu có 1.762.452 bản ghi và 17 biến. Quy mô này rất lớn, đủ mạnh cho phân tích mô tả chi tiết theo không gian–thời gian và điều kiện thời tiết.

Tên các biến

names(Accidents) 

##  [1] "ID"                "Source"            "Severity"         
##  [4] "Start_Time"        "End_Time"          "Start_Lat"        
##  [7] "Start_Lng"         "Distance(mi)"      "Street"           
## [10] "City"              "State"             "Weather_Timestamp"
## [13] "Temperature(F)"    "Humidity(%)"       "Pressure(in)"     
## [16] "Visibility(mi)"    "Weather_Condition"

Câu lệnh: (1) trả về một vector ký tự chứa tên tất cả các cột (biến) của bảng dữ liệu Accidents. Đây là cách nhanh để kiểm tra có bao nhiêu biến và đang đặt tên thế nào.

Kết quả: Kết quả hiển thị 17 tên biến, phản ánh đúng cấu trúc bộ US Accidents: ID, Source, Street, City, State, Start_Lat, Start_Lng, Start_Time, End_Time, Weather_Timestamp, Severity, Distance(mi), Temperature(F), Humidity(%), Pressure(in), Visibility(mi), Weather_Condition.

Kiểu dữ liệu từng biến

cls <- sapply(Accidents, \(x) paste(class(x), collapse="/"))
knitr::kable(matrix(sprintf("%s — %s", names(cls), cls), ncol=3, byrow=TRUE))

ID — character	Source — character	Severity — numeric
Start_Time — POSIXct/POSIXt	End_Time — POSIXct/POSIXt	Start_Lat — numeric
Start_Lng — numeric	Distance(mi) — numeric	Street — character
City — character	State — character	Weather_Timestamp — POSIXct/POSIXt
Temperature(F) — numeric	Humidity(%) — numeric	Pressure(in) — numeric
Visibility(mi) — numeric	Weather_Condition — character	ID — character

Câu lệnh: (1) Tạo vector chứa kiểu dữ liệu của từng cột trong bảng Accidents. (2) Ghép tên cột với kiểu dữ liệu thành chuỗi “tên – kiểu”, sắp xếp thành bảng 3 cột rồi hiển thị gọn bằng kable.

Kết quả: Chuỗi/ký tự (character): ID, Source, Street, City, State, Weather_Condition. Đây là các biến định danh/nhãn văn bản.

Số (numeric): Severity, Start_Lat, Start_Lng, Distance(mi), Temperature(F), Humidity(%), Pressure(in), Visibility(mi). Đây là các biến định lượng.

Ngày-giờ (POSIXct/POSIXt): Start_Time, End_Time, Weather_Timestamp. Kiểu này cho phép trích xuất giờ/thứ/tháng, sắp xếp theo thời gian, tính thời lượng, v.v.

Giới thiệu các biến và giải thích ý nghĩa từng biến

meaning <- c(ID = "Mã định danh duy nhất cho từng vụ tai nạn",
 Source = "Nguồn thu thập dữ liệu",
 Severity = "Mức độ nghiêm trọng (1 = nhẹ nhất, 4 = nghiêm trọng nhất)",
 Start_Time = "Thời điểm bắt đầu vụ tai nạn",
 End_Time = "Thời điểm kết thúc vụ tai nạn",
 Start_Lat = "Vĩ độ nơi bắt đầu vụ tai nạn",
 Start_Lng = "Kinh độ nơi bắt đầu vụ tai nạn",
 `Distance(mi)` = "Chiều dài đoạn đường bị ảnh hưởng (dặm)",
 Street = "Tên đường hoặc tuyến",
 City = "Thành phố nơi tai nạn ghi nhận",
 State = "Bang/khu vực",
 Weather_Timestamp = "Thời điểm thời tiết (timestamp quan trắc/thời điểm sự  kiện)",
 `Temperature(F)` = "Nhiệt độ (°F) tại thời điểm tai nạn",
 `Humidity(%)` = "Độ ẩm (%) tại thời điểm tai nạn",
 `Pressure(in)` = "Áp suất khí quyển (inch Hg)",
 `Visibility(mi)` = "Tầm nhìn (mile)",
 Weather_Condition = "Điều kiện thời tiết quan sát")
variable_meaning <- enframe(meaning, name = "Variable", value = "Meaning")
data_types <- map_chr(Accidents, ~ paste(class(.x), collapse = "|")) |>
  enframe(name = "Variable", value = "DataType")
variable_summary <- left_join(variable_meaning, data_types, by = "Variable")
knitr::kable(variable_summary, caption = "Tên biến, ý nghĩa và kiểu dữ liệu")

Tên biến, ý nghĩa và kiểu dữ liệu

Variable	Meaning	DataType
ID	Mã định danh duy nhất cho từng vụ tai nạn	character
Source	Nguồn thu thập dữ liệu	character
Severity	Mức độ nghiêm trọng (1 = nhẹ nhất, 4 = nghiêm trọng nhất)	numeric
Start_Time	Thời điểm bắt đầu vụ tai nạn	POSIXct|POSIXt
End_Time	Thời điểm kết thúc vụ tai nạn	POSIXct|POSIXt
Start_Lat	Vĩ độ nơi bắt đầu vụ tai nạn	numeric
Start_Lng	Kinh độ nơi bắt đầu vụ tai nạn	numeric
Distance(mi)	Chiều dài đoạn đường bị ảnh hưởng (dặm)	numeric
Street	Tên đường hoặc tuyến	character
City	Thành phố nơi tai nạn ghi nhận	character
State	Bang/khu vực	character
Weather_Timestamp	Thời điểm thời tiết (timestamp quan trắc/thời điểm sự kiện)	POSIXct|POSIXt
Temperature(F)	Nhiệt độ (°F) tại thời điểm tai nạn	numeric
Humidity(%)	Độ ẩm (%) tại thời điểm tai nạn	numeric
Pressure(in)	Áp suất khí quyển (inch Hg)	numeric
Visibility(mi)	Tầm nhìn (mile)	numeric
Weather_Condition	Điều kiện thời tiết quan sát	character

Câu lệnh: (1) Đoạn mã trên tạo bảng tổng hợp gồm tên biến, ý nghĩa và kiểu dữ liệu trong bộ dữ liệu Accidents. (2) Phần đầu khai báo vector meaning, chứa cặp giá trị tên biến và ý nghĩa bằng tiếng Việt. (3) để chuyển vector này thành bảng variable_meaning. (4) lấy kiểu dữ liệu (DataType) của từng biến trong Accidents và lưu vào bảng data_types. (5) ghép hai bảng lại theo cột Variable để tạo bảng variable_summary. (6) hiển thị bảng kết quả với tiêu đề “Tên biến, ý nghĩa và kiểu dữ liệu”, giúp mô tả cấu trúc dữ liệu một cách rõ ràng và dễ đọc.

Kết quả: Bảng trên trình bày 17 biến trong bộ dữ liệu cùng ý nghĩa và kiểu dữ liệu tương ứng.

ID: Mã định danh duy nhất cho mỗi vụ tai nạn, dạng ký tự (character), giúp tham chiếu và phát hiện trùng lặp.

Source: Nguồn thu thập dữ liệu (API giao thông, cảm biến…), dạng ký tự, không có đơn vị, dùng để so sánh chất lượng giữa các nguồn.

Severity: Mức độ nghiêm trọng từ 1–4, dạng số (numeric, ordinal), phản ánh thang đo mức độ ảnh hưởng của tai nạn.

Start_Time / End_Time: Thời điểm bắt đầu và kết thúc vụ tai nạn, dạng POSIXct, dùng để phân tích chu kỳ, kiểm tra khoảng thời gian ảnh hưởng.

Start_Lat / Start_Lng: Tọa độ vị trí tai nạn (vĩ độ, kinh độ), đơn vị độ thập phân, phục vụ cho phân tích không gian.

Distance(mi): Chiều dài đoạn đường bị ảnh hưởng, đơn vị mile, là biến định lượng phản ánh mức độ tác động của tai nạn.

Street / City / State: Tên đường, thành phố và bang nơi xảy ra tai nạn, dạng ký tự, phục vụ phân tích địa lý và bản đồ rủi ro.

Weather_Timestamp: Thời điểm quan trắc thời tiết, dạng POSIXct, dùng để đối chiếu sự trùng khớp giữa thời tiết và thời điểm tai nạn.

Temperature(F): Nhiệt độ tại thời điểm tai nạn, đơn vị °F, biến định lượng liên tục, có thể quy đổi sang °C = (°F − 32) × 5/9.

Humidity(%): Độ ẩm không khí tại hiện trường, đơn vị %, giá trị hợp lệ trong khoảng 0–100, ngoài phạm vi có thể là lỗi đo.

Pressure(in): Áp suất khí quyển, đơn vị inch thủy ngân (inHg), có thể quy đổi sang hPa, giá trị bất thường phản ánh lỗi cảm biến.

Visibility(mi): Tầm nhìn xa, đơn vị mile, thể hiện điều kiện quan sát của người lái, tầm nhìn thấp thường đi kèm mưa hoặc sương mù.

Weather_Condition: Mô tả điều kiện thời tiết định tính (Rain, Snow, Fog, Thunderstorm…), dạng ký tự, thường được gom nhóm để mô hình hóa.

Từ bảng trên có thể phân loại thành ba nhóm biến chính:

Biến định lượng: gồm Severity, Distance(mi), Temperature(F), Humidity(%), Pressure(in) và Visibility(mi) là các biến đo lường có giá trị số, phục vụ phân tích thống kê và đánh giá mức độ ảnh hưởng của tai nạn.

Biến định tính: gồm ID, Source, Street, City, State và Weather_Condition là các biến mô tả thông tin phân loại, nguồn dữ liệu và đặc điểm địa lý hoặc thời tiết.

Biến thời gian: gồm Start_Time, End_Time và Weather_Timestamp phản ánh thời điểm xảy ra sự kiện, phục vụ phân tích xu hướng, chu kỳ và dự báo tai nạn giao thông.

Kiểm tra chất lượng dữ liệu

Kiểm tra giá bị trị bị trùng lặp

Trong bộ US Accidents, mỗi dòng tương ứng một vụ tai nạn. Nếu có dòng trùng lặp, các thống kê về tần suất, mức độ (Severity) hay điều kiện thời tiết sẽ bị sai lệch. Do đó, kiểm tra trùng lặp là bước cần thiết để đảm bảo chất lượng và độ chính xác của dữ liệu trước khi phân tích.

sum(duplicated(Accidents))

## [1] 0

Câu lệnh: (1) sẽ kiểm tra các dòng trùng lặp trong toàn bộ bảng dữ liệu Accidents.

Kết quả: Kết quả trả về 0, nghĩa là không có hàng trùng hoàn toàn trong bảng mỗi dòng là một vụ tai nạn duy nhất theo toàn bộ 17 biến. Đây là tín hiệu tốt cho độ tin cậy của phần mô tả.

Kiểm tra phạm vi thời gian ghi nhận

Kiểm tra các mốc sớm–muộn nhất được ghi, xem có phủ đủ năm 2022 không, và phát hiện giá trị biên bất thường (ví dụ thời gian kết thúc vượt quá năm nghiên cứu).

range(Accidents$Start_Time, na.rm = TRUE)

## [1] "2022-01-01 00:02:00 UTC" "2022-12-31 23:59:03 UTC"

range(Accidents$End_Time,   na.rm = TRUE)

## [1] "2022-01-01 00:31:30 UTC" "2023-03-31 23:59:00 UTC"

Câu lệnh: (1) Trả về mốc sớm nhất–muộn nhất của Start_Time (bỏ qua NA) để biết khoảng thời gian phủ dữ liệu bắt đầu sự kiện. (2) Trả về mốc sớm nhất–muộn nhất của End_Time (bỏ qua NA) để biết khoảng thời gian phủ dữ liệu kết thúc sự kiện.

Kết quả: Start_Time trải từ 2022-01-01 00:02:00 UTC đến 2022-12-31 23:59:03 UTC ⇒ dữ liệu khởi phát sự kiện phủ trọn năm 2022; End_Time trải từ 2022-01-01 00:31:30 UTC đến 2023-03-31 23:59:00 UTC ⇒ có một số vụ kết thúc sau 2022, có thể do sự cố kéo dài/ghi nhận muộn.

Xử lý dữ liệu thô và mã hóa dữ liệu

Loại bỏ biến không cần thiết

Biến Weather_Timestamp ghi nhận thời điểm quan trắc thời tiết gần vụ tai nạn, nhưng thông tin này trùng với Start_Time và có nhiều giá trị thiếu. Các đặc trưng thời gian cần thiết đã được tách trực tiếp từ Start_Time, nên loại bỏ Weather_Timestamp giúp dữ liệu gọn hơn và tránh trùng lặp.

Accidents <- Accidents %>% dplyr::select(-Weather_Timestamp)

Câu lệnh: (1) dùng dplyr để cập nhật lại bảng dữ liệu, loại bỏ hoàn toàn cột Weather_Timestamp và giữ nguyên các cột còn lại.

Kiểm tra & mô tả thiếu dữ liệu

Kiểm tra giá trị bị thiếu

Trước khi phân tích, cần kiểm tra dữ liệu thiếu (NA) để đảm bảo chất lượng và độ tin cậy. Nếu không phát hiện cột thiếu nhiều, các phép tính có thể tự loại bớt hàng, làm giảm số mẫu hoặc gây sai lệch kết quả, nhất là khi phần thiếu xuất hiện trong trường hợp đặc biệt như mưa lớn.

colSums(is.na(Accidents))

##                ID            Source          Severity 
##                 0                 0                 0 
##        Start_Time          End_Time         Start_Lat 
##                 0                 0                 0 
##         Start_Lng      Distance(mi)            Street 
##                 0                 0              7236 
##              City             State    Temperature(F) 
##                65                 0             38718 
##       Humidity(%)      Pressure(in)    Visibility(mi) 
##             41107             33134             41953 
## Weather_Condition 
##             38688

Câu lệnh: (1) Câu lệnh tính tổng số giá trị NA trong từng cột của bộ dữ liệu Accidents, giúp xác định biến nào có dữ liệu bị thiếu.

Kết quả: Kết quả kiểm tra giá trị khuyết cho thấy bộ dữ liệu có một số biến chứa nhiều giá trị bị thiếu (NA).

Biến City có một lượng nhỏ giá trị thiếu (65). Tuy nhiên, các biến Street, Temperature(F), Humidity(%), Pressure(in), Visibility(mi), Weather_Condition có số lượng NA rất lớn lần lượt là 7.236, 29.595, 38.718, 41.107, 33.134, 41.953 và 38.688 giá trị.

ID, Source, Severity, Start_Time, End_Time, Start_Lat, Start_Lng, State không có giá trị thiếu, cho thấy dữ liệu về thời gian và vị trí tai nạn được ghi nhận đầy đủ.

Thiếu dữ liệu ở City/Street thường do tai nạn xảy ra trên cao tốc hoặc khu vực không có tên đường. Nhóm biến thời tiết thiếu vì được ghép từ trạm khí tượng gần nhất, khi xa trạm, lệch thời điểm hoặc trạm gián đoạn đo, hệ thống sẽ để trống. Điều này cho thấy dữ liệu có thể bị thiếu cục bộ do sự cố cảm biến hoặc kết nối.

Tính tỷ lệ giá trị bị thiếu NA

colMeans(is.na(Accidents)) * 100

##                ID            Source          Severity 
##       0,000000000       0,000000000       0,000000000 
##        Start_Time          End_Time         Start_Lat 
##       0,000000000       0,000000000       0,000000000 
##         Start_Lng      Distance(mi)            Street 
##       0,000000000       0,000000000       0,410564373 
##              City             State    Temperature(F) 
##       0,003688044       0,000000000       2,196825786 
##       Humidity(%)      Pressure(in)    Visibility(mi) 
##       2,332375577       1,879994462       2,380376884 
## Weather_Condition 
##       2,195123612

Câu lệnh: (1) dùng để tính tỷ lệ phần trăm giá trị thiếu (NA) cho từng cột trong bộ dữ liệu Accidents bằng cách lấy trung bình số NA trên tổng số dòng rồi nhân 100.

Kết quả: Kết quả trên thể hiện tỷ lệ phần trăm giá trị thiếu (NA) của từng biến trong bộ dữ liệu:

Các biến ID, Source, Severity, Start_Time, End_Time, Start_Lat, Start_Lng, State không có giá trị thiếu , dữ liệu nền tảng đầy đủ. Còn các biến bị thiếu: City gần như đầy đủ, chỉ thiếu 0,0037%, Street thiếu khoảng 0,41%, Temperature(F) 2,20%, Humidity(%) 2,33%, Pressure(in) 1,88%, Visibility(mi) 2,38%, Weather_Condition 2,20%. ⇒ Nhìn chung, dữ liệu tai nạn được ghi nhận tốt, chỉ nhóm biến thời tiết cần xử lý hoặc xem xét bổ sung trước khi mô hình hóa liệu.

Vẽ biểu đồ thể hiện tỷ lệ thiếu NA

Vẽ biểu đồ tỷ lệ thiếu giúp nhìn nhanh biến nào thiếu nhiều/ít để ra quyết định xử lý (giữ, điền, hay loại bỏ). So với nhìn số tuyệt đối, phần trăm cho phép so sánh công bằng giữa các biến và dễ thấy nhóm thời tiết thường thiếu hơn nhóm thời gian–không gian.

miss_df <- enframe(colMeans(is.na(Accidents)) * 100,
                   name = "variable", value = "pct_NA")
ggplot(miss_df, aes(x = 1,
           y = fct_rev(fct_reorder(variable, pct_NA)),
           fill = pct_NA)) +
  geom_tile(width = .9, height = .85) +
  geom_text(aes(label = sprintf("%.1f%%", pct_NA)), size = 3.2) +
  scale_fill_gradient(low = "#e8f5e9", high = "#ef5350", name = "Tỷ lệ NA") +
  labs(title = "TỶ LỆ THIẾU NA THEO BIẾN", x = NULL, y = "Biến") +
  theme_minimal(base_size = 14) +
  theme(panel.grid = element_blank(),
        axis.text.x = element_blank(),
        axis.ticks = element_blank())

Câu lệnh: (1) Tính tỷ lệ thiếu NA theo phần trăm và chuyển thành bảng miss_df. (2) Dùng ggplot() vẽ biểu đồ ô màu (heatmap) thể hiện tỷ lệ thiếu của từng biến. (3) geom_tile() tạo ô màu, geom_text() hiển thị giá trị phần trăm trên ô. (4) Thang màu từ xanh nhạt → đỏ đậm biểu thị tỷ lệ thiếu tăng dần. (5) Tiêu đề và theme được tinh giản để nhấn mạnh nội dung chính.

Kết quả: Kết quả biểu đồ cho thấy phần lớn các biến trong bộ dữ liệu Accidents có tỷ lệ thiếu 0%, bao gồm Distance(mi), End_Time, ID, Severity, Source, Start_Lat, Start_Lng, Start_Time, State và City.

Một số biến có tỷ lệ thiếu nhẹ, trong khoảng 0,4% – 2,4%, chủ yếu thuộc nhóm thời tiết và địa lý chi tiết như Street (0,4%), Pressure(in) (1,9%), Weather_Condition (2,2%), Temperature(F) (2,2%), Humidity(%) (2,3%) và Visibility(mi) (2,4%). Đây đều là những biến quan trọng trong việc lý giải mức độ nghiêm trọng của tai nạn, do đó không nên loại bỏ hoàn toàn mà cần xử lý cẩn thận.

Vẽ biểu đồ phân phối cho từng biến định lượng NA

num_vars <- c( "Severity","Start_Lat","Start_Lng", "Distance(mi)","Temperature(F)","Humidity(%)",
  "Pressure(in)","Visibility(mi)")
df_num <- Accidents[, num_vars] %>%
  tidyr::drop_na(dplyr::all_of(num_vars))
make_hist <- function(dat, var){
  ggplot(dat, aes(.data[[var]])) +
    geom_histogram(aes(y = after_stat(density)),
                   bins = 30, fill = "#6ea8fe", color = "white", alpha = .7) +
    geom_density(color = "#e45756", linewidth = .8, na.rm = TRUE) +
    labs(x = var, y = NULL) +                    
    theme_minimal(base_size = 9) +
    theme(
      axis.text    = element_text(size = 8),
      axis.title.x = element_text(margin = margin(t = 2)),
      plot.margin  = margin(3, 6, 3, 6))
}
plots  <- lapply(num_vars, function(v) make_hist(df_num, v))
p_all  <- patchwork::wrap_plots(plots, ncol = 3) +
  patchwork::plot_annotation(
    title    = "PHÂN PHỐI CÁC BIẾN ĐỊNH LƯỢNG",
    subtitle = "Histogram phủ đường mật độ, đã bỏ NA",
    theme = theme(
      plot.title    = element_text(face = "bold", size = 16, hjust = .5, margin = margin(b = 6)),
      plot.subtitle = element_text(size = 11,  hjust = .5, margin = margin(b = 8)) )) &
  theme(plot.margin = margin(3, 6, 3, 6))         
p_all

Câu lệnh: (1) Chọn các biến định lượng và lọc thành bảng df_num. (2) Tạo hàm make_hist() để vẽ histogram kèm đường mật độ cho từng biến. (3) geom_histogram() hiển thị tần suất, geom_density() biểu diễn phân bố trơn. (4) Dùng lapply() lặp qua các biến và wrap_plots() ghép tất cả biểu đồ vào một khung. (5) Tiêu đề và phụ đề giúp nêu rõ rằng biểu đồ thể hiện phân phối của các biến định lượng sau khi loại bỏ NA.

Kết quả:

Start_Lat: đa đỉnh rõ rệt quanh các vĩ độ ~30–40°B, điều này phản ánh cụm dân cư/đường xá chứ không phải nhiễu. Dùng để mô tả theo bang/vùng hoặc làm bản đồ là hợp lý.

Start_Lng: cũng đa đỉnh, với ba cụm tương ứng bờ Đông (~−75 đến −85), Midwest (~−95) và bờ Tây (~−120 đến −115). Mẫu hình này xác nhận phạm vi địa lý của dữ liệu và gợi ý phân tích theo khu vực.

Distance(mi): lệch phải rất mạnh và có nhiều giá trị 0.

Temperature(F) và Pressure có dạng phân phối gần chuẩn, tương đối đối xứng quanh giá trị trung tâm. Điều này cho thấy dữ liệu của hai biến này ổn định, không xuất hiện giá trị ngoại lai lớn.

=> Vì vậy, việc sử dụng giá trị trung bình (Mean) để thay thế các giá trị bị thiếu là phù hợp, giúp duy trì đặc trưng của phân phối ban đầu.

Các biến Humidity(%) và Visibility(mi) đều có phân phối lệch trái, phải rõ rệt và xuất hiện nhiều giá trị ngoại lai.

=> Vì vậy, dùng trung bình có thể làm méo phân phối do ảnh hưởng của ngoại lai. Do đó, nhóm chọn trung vị (Median) để thay thế giá trị thiếu, vì trung vị phản ánh trung tâm thực tế và ít bị ảnh hưởng bởi các điểm cực trị.

Riêng biến Severity là biến định tính dạng rời rạc. Do không mang tính liên tục, biến này được xử lý bằng giá trị xuất hiện nhiều nhất (mode) để đảm bảo tính đại diện cho nhóm dữ liệu phổ biến nhất.

Kiểm tra độ lệch phân phối các biến NA

Định lượng hình dạng phân phối của các biến quan trọng để quyết định cách xử lý NA và chọn thống kê phù hợp. Hai chỉ số chuẩn là skewness (độ lệch trái/phải) và kurtosis (độ nặng đuôi). Nếu phân phối gần chuẩn (skew≈0, kurtosis≈0) thì có thể điền mean; nếu lệch mạnh/đuôi nặng thì nên ưu tiên median, biến đổi log, hoặc điền theo nhóm.

vars <- c("Temperature(F)","Humidity(%)","Pressure(in)","Visibility(mi)","Distance(mi)")
res <- t(sapply(Accidents[vars],
                \(x) c(skewness = skewness(x, na.rm = TRUE),
                      kurtosis = kurtosis(x, na.rm = TRUE))))
data.frame(variable = rownames(res), res, row.names = NULL)

Câu lệnh: (1) Chọn các biến cần đo (vars). (2) Với từng biến, tính skewness và kurtosis (bỏ NA), rồi transpose để mỗi biến là một dòng. (3) Gộp thành data.frame có cột variable, skewness, kurtosis để so sánh hình dạng phân phối và quyết định xử lý NA/biến đổi (mean vs. median/log).

Kết quả:

Temperature(F) có skewness ≈ −0.7, kurtosis ≈ 0.25, gần chuẩn nên có thể điền Mean.

Humidity(%) gần đối xứng nhưng bị chặn ở biên, thích hợp dùng Median.

Pressure(in) lệch mạnh, đuôi dày, nhạy với ngoại lai nên điền Median đáng tin cậy hơn Mean.

Visibility(mi) lệch trái và đuôi dày, nên dùng Median để tránh méo dữ liệu.

Distance(mi) lệch phải cực mạnh (skew ≈ 13), nên cũng điền Median.

→ Tóm lại, nhóm dùng Mean cho biến gần chuẩn và Median cho các biến lệch hoặc có ngoại lai.

Xử lý thiếu dữ liệu

Điền giá trị mean cho Temperature

Biến Temperature(F) được lựa chọn điền khuyết bằng giá trị trung bình: bảo toàn kỳ vọng của biến, giữ nguyên cỡ mẫu cho các thống kê mô tả và trực quan hóa theo thời gian–không gian.

Accidents$`Temperature(F)`[is.na(Accidents$`Temperature(F)`)] <-mean(Accidents$`Temperature(F)`, na.rm = TRUE)

Câu lệnh: (1) Thay tất cả NA của biến Temperature(F) bằng giá trị trung bình của chính biến này (tính với na.rm = TRUE). → Sau khi chạy, cột không còn NA; cách điền này phù hợp vì Temperature(F) gần chuẩn.

Điền giá trị median cho các biến định lượng khác

Điền khuyết bằng trung vị cho mọi biến định lượng có NA, ngoại trừ Temperature(F) biến đã được chọn điền trung bình. Việc dùng median cho các biến số còn lại là hợp lý với bộ dữ liệu này vì nhiều phân phối lệch/đuôi nặng hoặc bị trần.

num_vars <- names(Accidents)[sapply(Accidents, is.numeric) & !(names(Accidents) %in% c("Temperature(F)"))]
for (v in num_vars) { Accidents[[v]][is.na(Accidents[[v]])] <- median(Accidents[[v]], na.rm = TRUE)}

Câu lệnh: (1) Chọn danh sách tên các biến số (numeric) trong Accidents và loại trừ biến Temperature(F). (2) Với mỗi biến còn lại, thay tất cả NA bằng trung vị (Median) của chính biến đó (na.rm = TRUE).

Điền giá trị mode cho biến định tính

Các biến định tính dạng chuỗi (như Street, City, Weather_Condition, …) có tỉ lệ thiếu rất nhỏ. Với biến danh mục, chiến lược đơn giản và nhất quán là điền bằng giá trị xuất hiện nhiều nhất (mode): vừa giữ nguyên cỡ mẫu, vừa không giả định khoảng cách giữa các mức như khi dùng số.

cat_vars <- names(Accidents)[sapply(Accidents, is.character)]
for (v in cat_vars) {mode_value <- names(sort(table(Accidents[[v]]), decreasing = TRUE))[1]
Accidents[[v]][is.na(Accidents[[v]])] <- mode_value}

Câu lệnh: (1) Chọn danh sách các biến dạng ký tự trong dữ liệu để xử lý. (2) Với từng biến, tính mode bằng cách lấy giá trị xuất hiện nhiều nhất. (3) Thay mọi NA của biến đó bằng mode tương ứng.

Kiểm tra lại xem còn giá trị NA

colSums(is.na(Accidents))

##                ID            Source          Severity 
##                 0                 0                 0 
##        Start_Time          End_Time         Start_Lat 
##                 0                 0                 0 
##         Start_Lng      Distance(mi)            Street 
##                 0                 0                 0 
##              City             State    Temperature(F) 
##                 0                 0                 0 
##       Humidity(%)      Pressure(in)    Visibility(mi) 
##                 0                 0                 0 
## Weather_Condition 
##                 0

Câu lệnh: (1) Lệnh trên đếm tổng số giá trị thiếu (NA) ở từng cột của bảng Accidents

Kết quả: Sau khi xử lý giá trị thiếu NA của các biến thì kết quả kiểm tra lại là 0.

Chuẩn hoá thời gian

Đinh dạng ngày giờ

Accidents <- Accidents |> dplyr::mutate(across(c(Start_Time, End_Time),
   ~ as.POSIXct(.x, format = "%Y-%m-%d %H:%M:%S", tz = "UTC")))

Câu lệnh: (1) Áp dụng mutate(across(…)) để xử lý đồng thời hai cột thời gian (Start_Time, End_Time). (2) Ép mỗi cột về POSIXct bằng mẫu “%Y-%m-%d %H:%M:%S” và múi giờ UTC.

Phân loại ngày tháng năm xảy ra tai nạn

Accidents <- Accidents %>% mutate(Date = as.Date(Start_Time))

Câu lệnh: (1) Lệnh trên tạo thêm cột Date bằng cách lấy phần ngày từ Start_Time thông qua as.Date().

Thêm cột giờ khi bắt đầu xảy ra vụ tai nạn

Mục tiêu là rút trích giờ trong ngày để phân tích nhịp tai nạn theo giờ (giờ cao điểm sáng/chiều, ban đêm…).

Accidents <- Accidents %>% mutate(Hour = as.numeric                       (format(Start_Time, "%H")))

Câu lệnh: (1) Tạo cột Hour là giờ trong ngày (0–23) được trích từ Start_Time và ép về số để dễ nhóm/tính.

Thêm cột ngày trong tuần xảy ra vụ tai nạn

Thêm ngày trong tuần của vụ tai nạn để phân tích nhịp tai nạn theo lịch làm việc/cuối tuần (so sánh weekday vs weekend, xem ngày nào cao nhất).

Accidents <- Accidents %>% mutate(Day = weekdays(Date))

Câu lệnh: (1) Tạo cột Day chứa tên thứ trong tuần từ biến Date bằng weekdays().

Gom nhóm & chuyển đổi biến

Thêm cột cột nhiệt độ mới theo độ °C

Chuẩn hoá về hệ SI để báo cáo/so sánh dễ hiểu (đa số độc giả quen °C).

Accidents <- Accidents %>% mutate(temp_c = (`Temperature(F)` - 32) * 5/9)

Công thức đổi: °C = (°F - 32) * 5/9.

Câu lệnh: (1) thêm cột temp_c và tính độ °C bằng công thức đổi.

Phân loại thời điểm xảy ra tai nạn trong ngày

Gộp 24 giờ thành vài “ca” (đêm–sáng–chiều–tối) giúp xem xu hướng dễ đọc hơn, tránh rời rạc theo từng giờ, và so sánh tần suất/tỷ lệ Severity theo thời điểm trong ngày một cách gọn, rõ.

Accidents <- Accidents %>%
  mutate(time_period = cut(Hour, breaks = c(0, 6, 12, 18, 24),                                               labels = c("Đêm", "Sáng", "Chiều","Tối"),                                   include.lowest = TRUE)) %>%
  mutate(time_period = factor(time_period,                                                         levels = c("Đêm","Sáng","Chiều","Tối"),                                     ordered = TRUE))

Câu lệnh: (1) Tạo biến time_period bằng cách chia biến Hour thành 4 khoảng có nhãn “Đêm – Sáng – Chiều – Tối”. (2) Dùng cut() để phân loại giờ theo nhóm, giúp rút gọn 24 mức thành 4 ca dễ so sánh. (3) Chuyển time_period thành factor có thứ tự, đảm bảo hiển thị đúng trình tự thời gian khi vẽ biểu đồ hoặc thống kê.

Phân loại điều kiện thời tiết lúc xảy ra tai nạn

Weather_Condition trong dữ liệu có rất nhiều mức chi tiết. Nếu giữ nguyên, bảng chéo/biểu đồ sẽ rối và khó diễn giải. Gộp chúng thành ít nhóm khí tượng có ý nghĩa giúp mô tả xu hướng mức độ nghiêm trọng theo điều kiện thời tiết rõ ràng hơn và ổn định thống kê hơn.

Accidents <- Accidents %>% mutate(
weather_simple = case_when(
str_detect(Weather_Condition, regex ("Rain|Thunderstorm", ignore_case=TRUE)) ~ "Mưa",
str_detect(Weather_Condition, regex ("Snow|Sleet|Ice", ignore_case=TRUE)) ~ "Tuyết",
str_detect(Weather_Condition, regex ("Fog|Mist|Haze", ignore_case=TRUE)) ~ "Sương mù",
str_detect(Weather_Condition, regex ("Cloud|Overcast", ignore_case=TRUE)) ~ "Nhiều mây",
str_detect(Weather_Condition, regex ("Clear|Fair", ignore_case=TRUE)) ~ "Trời quang", TRUE ~ "Khác"))

Câu lệnh: (1) Tạo biến weather_simple để gom nhóm điều kiện thời tiết. (2) Nếu Weather_Condition khớp regex Rain|Thunderstorm ⇒ gán “Mưa”. (3) Khớp Snow|Sleet|Ice ⇒ “Tuyết”. (4) Khớp Fog|Mist|Haze ⇒ “Sương mù”. (5) Khớp Cloud|Overcast ⇒ “Nhiều mây”. (6) Khớp Clear|Fair ⇒ “Trời quang”. (7) Mọi trường hợp còn lại ⇒ “Khác”.

Các thống kê cơ bản về bộ dữ liệu

Thống kê mô tả biến định tính

Xác định biến định tính

cat_vars <- names(Accidents)[sapply(Accidents, inherits, c("character","factor"))]
cat_vars

## [1] "ID"                "Source"            "Street"           
## [4] "City"              "State"             "Weather_Condition"
## [7] "Day"               "time_period"       "weather_simple"

Câu lệnh: (1) trả về TRUE nếu cột là chuỗi hoặc factor → đây là hai kiểu dữ liệu tiêu biểu của biến định tính. Kết quả là một vector logic; đặt trong chỉ số của names(Accidents)[ ... ] để lấy tên các cột thoả điều kiện. (2) cho ta danh sách biến định tính hiện có.

Thống kê giá trị xuất hiện nhiều nhất của mỗi dữ liệu biến định tính

freq_all <- Accidents %>%
  dplyr::select(where(~is.character(.x) || is.factor(.x)), -ID, -Source) %>%
  tidyr::pivot_longer(everything(), names_to="Bien", values_to="Gia_tri") %>%
  count(Bien, Gia_tri, name="n") %>%
  group_by(Bien) %>% mutate(Ty_le_pct = 100*n/sum(n)) %>% ungroup() %>%
  mutate(Ty_le_pct_VN = format(round(Ty_le_pct,2), decimal.mark=",", nsmall=2))
mode_tbl <- freq_all %>%
  group_by(Bien) %>% slice_max(n, n=1, with_ties=FALSE) %>% ungroup() %>%
  transmute(Bien,
            Gia_tri_pho_bien = Gia_tri,
            So_lan = format(n, big.mark=".", decimal.mark=",", scientific=FALSE))
mode_tbl

Câu lệnh: (1–6) Tạo bảng tần suất cho biến định tính: chọn các cột dạng chữ/factor (loại ID, Source) → chuyển dữ liệu sang dạng dài (cột Biến và Giá trị) → đếm số lần từng cặp xuất hiện → theo từng Biến tính tỉ lệ % và định dạng % kiểu VN. (7–11) Rút mode cho từng biến: từ bảng tần suất ở trên, nhóm theo Biến và chọn hàng có n lớn nhất (không cho phép hòa) → chỉ giữ 3 cột: Biến, Giá trị phổ biến nhất, Số lần (định dạng có chấm ngăn cách nghìn, không scientific). (12) Hiển thị kết quả.

Kết quả: Miami là thành phố có số vụ tai nạn nhiều nhất (64.609), phản ánh mật độ ghi nhận dữ liệu cao hơn chứ không hẳn là nơi nguy hiểm nhất. Thứ Sáu ghi nhận nhiều vụ nhất (~307.733), phù hợp với thời điểm lưu thông cao cuối tuần. Bang California (CA) dẫn đầu với 375.913 bản ghi, thể hiện dữ liệu tập trung tại các bang đông dân. Tuyến đường I-95 S xuất hiện nhiều nhất (23.999), gợi ý đây là điểm nóng giao thông. Về thời tiết, điều kiện “Trời quang” chiếm đa số (904.354), cho thấy phần lớn tai nạn xảy ra trong điều kiện bình thường chứ không do thời tiết xấu. Khung thời gian “Chiều” (758.195) cũng nổi bật, trùng với giờ tan tầm – thời điểm lưu lượng xe cao nhất.

Thống kê tần suất Top 10 thành phố xảy ra nhiều tai nạn nhất

top10_city <- Accidents %>%
  count(City, name = "n") %>%
  mutate(ty_le = round(n / sum(n) * 100, 2)) %>%
  arrange(desc(n)) %>% slice_head(n = 10) %>%
  mutate(n = format(n, big.mark = ".", scientific = FALSE),
    ty_le = gsub("\\.", ",", format(ty_le, nsmall = 2)))
top10_city

Câu lệnh: (1) Khởi tạo pipeline từ Accidents và đặt tên kết quả top10_city. (2) Đếm số bản ghi theo City, tạo cột n. (3) Tính ty_le = tỷ lệ % của từng thành phố (làm tròn 2 số). (4) Sắp xếp giảm dần theo n và lấy Top-10. (5) Định dạng n có dấu chấm ngăn cách nghìn. (6) Định dạng ty_le theo chuẩn VN (dấu phẩy, 2 chữ số thập phân). (7) Hiển thị bảng top10_city.

Kết quả: Miami dẫn đầu với 64.609 vụ (3,67%), theo sau là Orlando và Los Angeles. Các thành phố lớn, đông dân hoặc có hạ tầng giao thông dày đặc thường ghi nhận số vụ cao hơn, phản ánh mật độ lưu thông và mức độ ghi nhận dữ liệu chứ không hẳn mức độ nguy hiểm.

Thống kê tần suất Top 10 đường xảy ra nhiều tai nạn nhất

top10_street <- Accidents %>%
  count(Street, name = "n") %>%                        
  mutate(ty_le = round(n / sum(n) * 100, 2)) %>%    
  arrange(desc(n)) %>%                              
  slice_head(n = 10) %>%
  mutate( n = format(n, big.mark = ".", scientific = FALSE),
    ty_le = gsub("\\.", ",", format(ty_le, nsmall = 2)))  
top10_street

Câu lệnh: (1) Tạo bảng top10_street từ dữ liệu gốc. (2) Đếm số vụ theo Street → cột đếm tên n. (3) Tính ty_le = tỷ trọng % của mỗi bang trong tổng, làm tròn 2 chữ số. (4) Sắp xếp n giảm dần và lấy Top10_street. (5) Hiển thị 10 dòng cao nhất. (6) Định dạng n có dấu chấm ngăn cách nghìn, tắt scientific và (6) Định dạng ty_le giữ 2 chữ số thập phân và dùng dấu , kiểu Việt Nam.

Nhận xét: I-95 là tuyến có nhiều tai nạn nhất, gồm cả các nhánh I-95 South, I-95 North và I-95 chung, cho thấy trục giao thông dọc bờ Đông là điểm nóng tai nạn lớn nhất toàn mạng lưới. Các tuyến I-5 và I-10 cũng nằm trong nhóm đầu, phản ánh mật độ phương tiện cao trên các hành lang bờ Tây và vành đai Nam Hoa Kỳ. Tỷ lệ tai nạn dao động khoảng 0,46–1,36%, khá nhỏ do dữ liệu phân tán trên nhiều tuyến, nhưng vẫn thể hiện sự tập trung đáng kể ở một số cao tốc chính.

Thống kê tần suất Top 10 bang xảy ra nhiều tai nạn nhất

top10_bang <- Accidents %>%
  count(State, name = "n") %>%
  arrange(desc(n)) %>% slice_head(n = 10) %>%
  mutate(ty_le = round(100*n/sum(n), 2)) %>%
  mutate(
    n = format(n, big.mark = ".", scientific = FALSE),
    ty_le = gsub("\\.", ",", format(ty_le, nsmall = 2))) 
top10_bang

Câu lệnh: (1) Tạo bảng top10_bang từ dữ liệu gốc. (2) Đếm số vụ theo State. (3) Sắp xếp n giảm dần và lấy Top-10 bang. (4) Tính tỷ lệ, làm tròn 2 chữ số. (5) Chuẩn hoá hiển thị: (6) Định dạng n có dấu chấm ngăn cách nghìn, tắt scientific và (6) Định dạng ty_le giữ 2 chữ số thập phân và dùng dấu , kiểu Việt Nam.

Kết quả: CA (30,12%), FL (21,08%) chiếm >51% toàn bộ mẫu → phân bố rất lệch (concentration), hai bang đông dân/đô thị lớn kéo tổng số vụ lên mạnh. Nhóm trung vị gồm VA (7,96%), TX (7,65%), NY (7,59%), SC (6,80%), PA (6,48%), NC (5,78%): mức khá sát nhau → phản ánh quy mô giao thông lớn nhưng không cực đoan như CA/FL. Đuôi nhỏ: NJ (3,38%), AZ (3,16%) → đóng góp ít hơn nhiều, phù hợp dân số/diện tích hoặc mạng lưới đường thấp hơn trong mẫu.

Thống kê tần suất số vụ tai nạn theo nhóm điều kiện thời tiết

weather_summary <- Accidents %>%
  count(weather_simple, name = "So_vu") %>%
  arrange(desc(So_vu)) %>%
  mutate(Ty_le = round(So_vu / sum(So_vu) * 100, 2)) %>%
  mutate(
    So_vu = format(So_vu, big.mark = ".", scientific = FALSE),
    Ty_le = gsub("\\.", ",", format(Ty_le, nsmall = 2)))
weather_summary

Câu lệnh: (1) Tạo bảng weather_summary từ dữ liệu Accidents để tổng hợp số vụ theo điều kiện thời tiết. (2) đếm số vụ tai nạn theo từng nhóm điều kiện thời tiết đơn giản, lưu kết quả vào cột So_vu. (3) sắp xếp các nhóm thời tiết theo thứ tự giảm dần số vụ. (4) tính tỷ lệ phần trăm số vụ mỗi loại thời tiết so với tổng số vụ và làm tròn 2 chữ số. (5) bắt đầu một mutate() mới để định dạng hiển thị. (6) định dạng số vụ có dấu chấm ngăn nghìn, không dùng ký hiệu khoa học. (7) chuyển dấu thập phân từ “.” sang “,” cho đúng kiểu Việt Nam. (8) xuất kết quả.

Kết quả:

1. Trời quang chiếm tỷ lệ cao nhất (52,03%), nghĩa là hơn một nửa số vụ tai nạn xảy ra trong điều kiện thời tiết tốt, tầm nhìn xa rõ. Điều này cho thấy tai nạn không chỉ phụ thuộc vào yếu tố thời tiết, mà còn liên quan đến hành vi lái xe, mật độ phương tiện và hạ tầng giao thông.

2. Nhiều mây chiếm 35,58%, là điều kiện phổ biến thứ hai — vẫn là thời tiết tương đối ổn định nhưng có thể gây giảm nhẹ ánh sáng, ảnh hưởng tầm nhìn và phản ứng của người lái.

3. Các hiện tượng bất lợi như mưa (5,70%), tuyết (2,60%) và sương mù (2,06%) chỉ chiếm tỷ trọng nhỏ. Tuy nhiên, dù ít xảy ra, chúng thường làm tăng mức độ nghiêm trọng của tai nạn do trơn trượt, hạn chế tầm nhìn hoặc mất kiểm soát phương tiện.

4. Nhóm “Khác” chỉ 2,03% phản ánh các trạng thái thời tiết hiếm gặp, có thể bao gồm bão, gió mạnh hoặc điều kiện đặc biệt khác.

Thống kê tần suất ngày trong tuần xảy ra tai nạn

day_summary <- Accidents %>%
  count(Day, name = "So_vu") %>%                     
  mutate(Ty_le = round(So_vu / sum(So_vu) * 100, 2)) %>%  
  arrange(desc(So_vu)) %>%
  mutate(
    So_vu = format(So_vu, big.mark = ".", scientific = FALSE),
    Ty_le = gsub("\\.", ",", format(Ty_le, nsmall = 2)))                              
day_summary

Câu lệnh: (1) Tạo bảng day_summary từ dữ liệu Accidents để tổng hợp số vụ theo ngày. (2) đếm số vụ tai nạn theo từng nhóm day. (3) tính tỷ lệ phần trăm số vụ và làm tròn 2 chữ số. (4) sắp xếp các nhóm theo thứ tự giảm dần số vụ. (5) bắt đầu một mutate() mới để định dạng hiển thị. (6) định dạng số vụ có dấu chấm ngăn nghìn, không dùng ký hiệu khoa học. (7) chuyển dấu thập phân từ “.” sang “,” cho đúng kiểu Việt Nam. (8) xuất kết quả.

Kết quả: Biểu đồ thể hiện tỷ lệ tai nạn giao thông theo ngày trong tuần cho thấy mức độ chênh lệch rõ rệt giữa các ngày. Các ngày Thứ 4, 5 và 6 chiếm tỷ trọng cao nhất (khoảng 16–17%) tương ứng với giai đoạn cao điểm đi lại trong tuần làm việc. Ngược lại, Chủ nhật có tỷ lệ thấp nhất (~8,8%), có thể do lưu lượng giao thông giảm và ít xe thương mại lưu thông. Tổng thể, xu hướng cho thấy tai nạn tập trung vào giữa và cuối tuần, phản ánh mật độ giao thông và áp lực di chuyển tăng dần từ đầu tuần tới cuối tuần. Điều này gợi ý cơ quan quản lý nên tăng cường kiểm soát và tuyên truyền an toàn giao thông vào các ngày giữa tuần, đặc biệt là Thứ 5–Thứ 6, khi rủi ro xảy ra cao nhất.

Thống kê tần suất thời điểm trong ngày xảy ra tai nạn

time_summary <- Accidents %>%
  count(time_period, name = "So_vu") %>%
  mutate(Ty_le = round(So_vu / sum(So_vu) * 100, 2)) %>%
  arrange(desc(So_vu))%>%
  mutate(
    So_vu = format(So_vu, big.mark = ".", scientific = FALSE),
    Ty_le = gsub("\\.", ",", format(Ty_le, nsmall = 2)))                  
time_summary

Câu lệnh: (1) Tạo bảng time_summary từ dữ liệu Accidents để tổng hợp số vụ theo thời điểm trong ngày. (2) đếm số vụ tai nạn theo từng nhóm theo thời điểm trong ngày, lưu kết quả vào cột So_vu. (3) tính tỷ lệ phần trăm số vụ so với tổng số vụ và làm tròn 2 chữ số. (4) sắp xếp các nhóm theo thứ tự giảm dần số vụ. (5) bắt đầu một mutate() mới để định dạng hiển thị. (6) định dạng số vụ có dấu chấm ngăn nghìn, không dùng ký hiệu khoa học. (7) chuyển dấu thập phân từ “.” sang “,” cho đúng kiểu Việt Nam. (8) xuất kết quả. Kết quả:

1. Thứ Sáu có số vụ tai nạn cao nhất (17,46%), cho thấy đây là ngày rủi ro nhất trong tuần — thường gắn với lưu lượng giao thông tăng mạnh khi kết thúc tuần làm việc, người dân di chuyển nhiều hơn (đi chơi, về quê, vận tải hàng hóa…).

2. Từ Thứ Tư đến Thứ Năm (16,31–16,60%) cũng có tỷ lệ cao, thể hiện xu hướng gia tăng tai nạn vào giữa và cuối tuần làm việc, khi cường độ giao thông cao và mức độ tập trung của người lái có thể giảm dần.

3. Thứ Hai và Thứ Ba có tỷ lệ thấp hơn (14,00% và 15,55%), phù hợp với giai đoạn đầu tuần, khi nhịp độ giao thông ổn định hơn.

4. Cuối tuần (Thứ Bảy – 11,30%, Chủ Nhật – 8,79%) có ít tai nạn hơn rõ rệt, phản ánh việc giảm lưu lượng xe cộ và áp lực giao thông trong ngày nghỉ.

=> Tai nạn xảy ra nhiều nhất vào cuối tuần làm việc (đặc biệt là Thứ Sáu), giảm mạnh vào hai ngày cuối tuần, cho thấy yếu tố hành vi di chuyển và mật độ giao thông là nguyên nhân chủ đạo ảnh hưởng đến tần suất tai nạn.

Thống kê mô tả biến định lượng

Thống kê mô tả biến Severity

summary(Accidents$Severity)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   1,000   2,000   2,000   2,072   2,000   4,000

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “Severity” (mức độ nghiêm trọng của tai nạn).

Kết quả: Phần lớn các vụ tai nạn có mức độ nghiêm trọng = 2, thể hiện qua median, 1st Qu. và 3rd Qu. đều bằng 2.mGiá trị mean ≈ 2,07 cho thấy phân phối hơi lệch nhẹ về phía các mức nặng hơn. Mức độ nghiêm trọng cao nhất là 4, nhưng hiếm gặp. → Nhìn chung, đa số vụ tai nạn có mức độ trung bình, còn các vụ cực nhẹ (1) hoặc rất nặng (4) chỉ chiếm tỷ lệ nhỏ.

Thống kê mô tả biến Distance(mi)

summary(Accidents$`Distance(mi)`)

##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
##   0,0000   0,0570   0,2600   0,9264   0,9780 336,5700

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “Distance(mi)” (chiều dài đoạn đường bị ảnh hưởng).

Kết quả: Biến Distance(mi) có phân phối lệch phải rất mạnh: 50% vụ tai nạn ảnh hưởng không quá 0,26 mile và 75% không quá 0,98 mile, cho thấy phần lớn sự cố chỉ diễn ra trong phạm vi ngắn, mang tính cục bộ tại các giao lộ hoặc khu dân cư đông đúc. Giá trị trung bình 0,93 mile cao hơn trung vị nhiều lần, phản ánh sự tồn tại của một số ít vụ có phạm vi ảnh hưởng rất lớn kéo đuôi phải của phân phối. Ở đầu trên, Max = 336,57 mile vượt xa Q3 (tỷ lệ Max/Q3 ≈ 344) là dấu hiệu rõ rệt của ngoại lệ, có thể xuất phát từ các vụ phong tỏa đường dài, tai nạn dây chuyền hoặc lỗi ghi nhận dữ liệu.

Thống kê mô tả biến Temperature(C)

summary(Accidents$`temp_c`)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  -42,78   10,56   18,89   17,32   25,56   77,78

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “temp_c” (nhiệt độ theo °C).

Kết quả: Nhiệt độ tại thời điểm tai nạn dao động rộng từ –42,78°F đến 77,78°F, trung bình 17,32°F và trung vị 18,89°F, cho thấy phân phối khá cân đối quanh mức lạnh–mát. Khoảng tứ phân vị (10,56–25,56°F) tập trung chủ yếu trong vùng nhiệt độ thấp, phản ánh nhiều vụ tai nạn xảy ra trong điều kiện thời tiết lạnh, có thể làm giảm tầm nhìn và độ bám đường.

Thống kê mô tả biến Humidity(%)

summary(Accidents$`Humidity(%)`)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1,00   45,00   63,00   62,06   81,00  100,00

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “Humidity(%)” (độ ẩm).

Kết quả: Độ ẩm dao động từ 1–100%, với trung bình 62,06% và trung vị 63%, cho thấy phân phối gần đối xứng và tập trung quanh mức ẩm cao. Khoảng 45–81% là dải phổ biến khi xảy ra 50% vụ tai nạn, phản ánh phần lớn sự cố diễn ra trong điều kiện ẩm ướt. Giá trị cực đại 100% thường gắn với sương mù hoặc mưa, trong khi mức 1% có thể là ngoại lệ hoặc điều kiện cực khô hiếm gặp.

Thống kê mô tả biến Pressure(in)

summary(Accidents$`Pressure(in)`)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0,00   29,21   29,72   29,38   29,96   58,63

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “Pressure(in)” (áp suất).

Kết quả: Áp suất tập trung chặt quanh mức chuẩn: Q1=29.21, Median=29.72, Mean=29.38, Q3=29.96 inHg—phù hợp dải điển hình gần mực biển (~29.92 inHg), cho thấy phân tán nhỏ và dữ liệu nhìn chung hợp lý. Hai cực trị Min=0.00 và Max=58.63 inHg là phi vật lý (ngoài khoảng thực tế ~28–31), nhiều khả năng lỗi đo/mã hóa (0 thay cho thiếu, 58.63 do nhập sai). Khi mô hình hóa, nên gán NA cho giá trị <25 hoặc >32 inHg (hoặc winsorize về rìa hợp lý).

Thống kê mô tả biến Visibility(mi)

summary(Accidents$`Visibility(mi)`)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   0,000  10,000  10,000   9,178  10,000 100,000

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “Visibility(mi)” (tầm nhìn).

Kết quả: Tầm nhìn (Visibility) có trung bình 9,178 dặm, trung vị 10 dặm, dao động từ 0 đến 100 dặm. Giá trị 0 và 100 là bất thường (0 có thể là sương mù dày hoặc lỗi ghi; 100 gần chắc là sai đơn vị/nhập liệu vì thang đo thực tế thường tối đa 10), trong đó phần lớn giá trị tập trung quanh 10 dặm – cho thấy đa số tai nạn xảy ra khi tầm nhìn khá tốt.

Thống kê mô tả biến Start_Lng

summary(Accidents$Start_Lng)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -124,54 -112,30  -84,50  -93,09  -80,09  -68,28

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “Start_Lng” (kinh độ).

Kết quả: Biến này thể hiện kinh độ (Start_Lng) của vị trí tai nạn, dao động từ -124,54 đến -68,28, tương ứng phạm vi toàn lãnh thổ Mỹ (từ bờ Tây – California đến bờ Đông – Maine). Trung vị -84,50 và trung bình -93,09 cho thấy phần lớn vụ tai nạn tập trung về phía Đông và Trung Mỹ – nơi có mật độ dân cư, đô thị và hệ thống giao thông cao hơn. Điều này hàm ý hoạt động kinh tế và di chuyển dày đặc hơn kéo theo rủi ro tai nạn cao hơn, đặc biệt tại các vùng công nghiệp và đô thị ven biển Đông.

Thống kê mô tả biến Hour

summary(Accidents$Hour)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0,00    8,00   14,00   12,73   17,00   23,00

Câu lệnh: (1) trong R dùng để tóm tắt mô tả biến “Hour”.

Kết quả: Biến Hour dao động từ 0–23h, với trung vị 14h và trung bình ≈12,7h, cho thấy tai nạn tập trung chủ yếu ban ngày (8h–17h) — trùng với khung giờ cao điểm đi làm, tan tầm. Điều này phản ánh mật độ giao thông và hoạt động kinh tế là yếu tố chính làm tăng rủi ro tai nạn, hơn là yếu tố thời tiết hay ban đêm.

Thống kê tần số một vài biến theo nhóm

Thời tiết và ca trong ngày

freq2_row <- Accidents %>%
  mutate(time_period = factor(time_period, levels = c("Sáng","Chiều","Tối","Đêm"))) %>%
  count(weather_simple, time_period, name = "n") %>%
  group_by(weather_simple) %>% mutate(p = n/sum(n)) %>% ungroup() %>%
  transmute(weather_simple, time_period,
            `n (%)` = sprintf("%s (%.2f%%)",
                              format(n, big.mark=".", decimal.mark=",", scientific=FALSE),
                              100*p)) %>%
  tidyr::pivot_wider(names_from = time_period, values_from = `n (%)`,
                     values_fill = "0 (0,00%)")
freq2_row

Câu lệnh: (1–3) Chuẩn hoá thứ tự time_period theo Sáng → Chiều → Tối → Đêm, rồi đếm số vụ theo từng cặp (weather_simple, time_period) để có tần suất n. (4) Tính tỷ trọng theo hàng trong mỗi trạng thái thời tiết. (5–7) Tạo cột hiển thị “n (xx,xx%)”: định dạng n có dấu chấm ngăn nghìn, phần trăm 2 chữ số thập phân dùng dấu phẩy. (8–10) Xoay bảng sang dạng rộng: mỗi cột là một buổi (Sáng/Chiều/Tối/Đêm), ô là chuỗi “n (xx,xx%)”; thiếu dữ liệu thì điền “0 (0,00%)”. (11) Xuất bảng kết quả tổng hợp.

Kết quả: Bảng cho thấy phần lớn tai nạn xảy ra vào buổi chiều (41–56%) và buổi sáng (29–45%), tập trung ở các điều kiện thời tiết trời quang và nhiều mây — hai nhóm chiếm hơn 70% tổng số vụ. Điều này phản ánh rằng tai nạn chủ yếu xảy ra trong giờ cao điểm giao thông, khi mật độ phương tiện lớn chứ không hoàn toàn do thời tiết xấu. → Về kinh tế, điều này hàm ý thiệt hại do tắc nghẽn và va chạm nhẹ trong đô thị là đáng kể, cần các biện pháp điều tiết giao thông theo khung giờ cao điểm và tăng cường cảnh báo an toàn ở điều kiện thời tiết bình thường.

Ca trong ngày và ngày trong tuần

freq_colpct <- Accidents %>%
  count(time_period, Day, name = "n") %>%
  group_by(time_period) %>% mutate(p = n/sum(n)) %>% ungroup() %>%
  transmute(time_period, Day,
            cell = sprintf("%s (%.2f%%)",
                           format(n, big.mark=".", decimal.mark=",", scientific=FALSE),
                           100*p)) %>%
  tidyr::pivot_wider(names_from = Day, values_from = cell,
                     values_fill = "0 (0,00%)")
freq_colpct

Câu lệnh: (1–2) Đếm số vụ tai nạn theo từng tổ hợp giữa ca trong ngày và ngày trong tuần. (3) Tính tỷ lệ phần trăm theo từng ca trong ngày: trong mỗi nhóm time_period, chia n cho tổng số của ca đó để ra tỷ trọng p, rồi bỏ nhóm. (4–6) Tạo cột mới hiển thị “n (xx,xx%)” với định dạng đẹp: có dấu chấm ngăn cách nghìn, dấu phẩy cho phần thập phân, phần trăm làm tròn 2 chữ số. (7–9) Xoay bảng sang dạng rộng, trong đó mỗi cột là một ngày trong tuần (Monday → Sunday), còn giá trị trong ô là chuỗi “n (xx,xx%)”; chỗ không có dữ liệu được điền “0 (0,00%)”.(10) Xuất bảng kết quả tổng hợp freq_colpct.

Kết quả: Kết quả cho thấy tai nạn tập trung cao nhất vào buổi chiều (≈19%) và tối (≈18%), thấp hơn vào buổi sáng (≈15%) và đêm (≈13–14%). Điều này phù hợp với thực tế khi lưu lượng xe tăng mạnh giờ tan tầm, tâm lý mệt mỏi và ánh sáng giảm làm tăng rủi ro va chạm. Cuối tuần (thứ 6, 7) vẫn duy trì tỷ lệ cao do hoạt động đi lại, giải trí và du lịch. → Về kinh tế, điều này phản ánh rủi ro chi phí xã hội gia tăng: tắc nghẽn giao thông, hư hại tài sản, gián đoạn lao động và logistics đô thị. Do đó, cần phân bổ nguồn lực giao thông hợp lý, tăng cường chiếu sáng, kiểm soát tốc độ và tuyên truyền an toàn vào các khung giờ cao điểm để giảm thiệt hại kinh tế – xã hội.

Ca trong ngày và mức độ nghiêm trọng

tab_sev3_shift <- Accidents %>%
  mutate(sev3 = ifelse(Severity >= 3, "Nặng (≥3)", "Nhẹ (<3)")) %>%
  count(time_period, sev3, name = "n") %>%
  group_by(time_period) %>% mutate(p = n/sum(n)) %>% ungroup() %>%
  transmute(time_period, sev3,
            cell = paste0(format(n, big.mark=".", decimal.mark=",", scientific=FALSE),
                          " (", format(round(100*p,2), decimal.mark=",", nsmall=2), "%)")) %>%
  tidyr::pivot_wider(names_from = sev3, values_from = cell, values_fill = "0 (0,00%)")
tab_sev3_shift

Câu lệnh: (1–2) Tạo biến sev3 hai mức (Nặng ≥3 / Nhẹ <3) từ Severity. (3) Đếm số vụ theo từng cặp ca trong ngày × mức sev3 để có tần suất n. (4) Trong mỗi ca, tính tỷ trọng. (5–7) Ghép chuỗi hiển thị “n (xx,xx%)”: n có dấu chấm ngăn nghìn; phần trăm làm tròn 2 số lẻ, dùng dấu phẩy. (8) Xoay bảng sang dạng rộng: mỗi cột là một mức sev3 (Nặng/ Nhẹ), mỗi ô là chuỗi “n (xx,xx%)”; thiếu dữ liệu thì điền “0 (0,00%)”. (9) Xuất bảng kết quả tổng hợp.

Kết quả: Kết quả cho thấy tai nạn nhẹ (<3) chiếm ưu thế tuyệt đối trong mọi khung giờ, khoảng 92–94%, trong khi tai nạn nặng (≥3) chỉ dao động 6–8%. Tỷ lệ tai nạn nặng cao nhất vào ban đêm (8,09%), tiếp theo là buổi tối (7,57%), sáng (7,21%), và thấp nhất buổi chiều (6,05%). Điều này phản ánh đặc trưng thực tế của giao thông: ban đêm và tối là thời điểm nguy hiểm hơn do tầm nhìn kém, mệt mỏi sau ngày làm việc và vi phạm tốc độ.

=> Cần ưu tiên can thiệp ban Đêm/Tối khi rủi ro tai nạn nặng cao nhất: tăng chiếu sáng, tuần tra tốc độ và kiểm soát rượu bia. Y tế khẩn cấp nên bố trí sẵn tại điểm nóng để giảm tử vong. Truyền thông cần nhấn mạnh “lái xe ban đêm rủi ro cao”, khuyến khích nghỉ ngơi và kiểm tra đèn – kính trước khi đi.

Ca trong ngày và nhiệt độ

tab_temp_band <- Accidents %>%
  mutate(temp_band = cut(temp_c,
                         breaks = c(-Inf, 10, 25, Inf),
                         labels = c("Thấp (≤10°C)", "Vừa (10–25°C)", "Cao (>25°C)"))) %>%
  count(temp_band, time_period, name = "n") %>%
  group_by(temp_band) %>% mutate(p = n/sum(n)) %>% ungroup() %>%
  transmute(temp_band, time_period,
            cell = paste0(format(n, big.mark=".", decimal.mark=",", scientific=FALSE),
                          " (", format(round(100*p,2), decimal.mark=",", nsmall=2), "%)")) %>%
  tidyr::pivot_wider(names_from = time_period, values_from = cell, values_fill = "0 (0,00%)")
tab_temp_band 

Câu lệnh: (1–2) Chia nhiệt độ temp_c thành 3 nhóm temp_band: Thấp (≤10°C), Vừa (10–25°C), Cao (>25°C) theo các ngưỡng đã đặt. (5) Đếm số vụ cho từng cặp (nhóm nhiệt độ × ca trong ngày) để thu được tần suất n. (6) Trong mỗi nhóm nhiệt độ, tính tỷ trọng. (7–9) Tạo cột hiển thị “n (xx,xx%)”: n định dạng có dấu chấm ngăn nghìn; phần trăm làm tròn 2 số lẻ và dùng dấu phẩy. (10) Xoay bảng về dạng rộng: mỗi cột là ca trong ngày, ô là chuỗi “n (xx,xx%)”; ô trống điền “0 (0,00%)”. (11) Xuất bảng kết quả tổng hợp.

Kết quả:

Cao (>25°C): Tập trung mạnh vào Chiều (63,84%), tiếp đến Sáng (24,40%); Đêm và Tối rất thấp (2,81% và 8,94%). → Thời tiết nóng chủ yếu xuất hiện buổi chiều.

Vừa (10–25°C): Phân bổ lệch về Chiều (38,99%), Sáng ~31,39%, Đêm/Tối ~15%. → Dải nhiệt độ dễ chịu vẫn nghiêng về ban ngày, đặc biệt buổi chiều.

Thấp (≤10°C): Sáng chiếm cao nhất (34,43%), tiếp theo Chiều (29,18%), Đêm (21,44%), Tối (14,94%). → Nhiệt thấp xuất hiện nhiều hơn vào sáng sớm/đêm, phù hợp đặc tính lạnh về đêm và buổi sáng.

Bất kể dải nhiệt, buổi chiều chiếm tỷ trọng lớn nhất (đặc biệt khi nhiệt cao); đêm/tối giảm mạnh khi nhiệt tăng. Điều này gợi ý rủi ro liên quan nhiệt/ánh nắng tập trung buổi chiều, còn lạnh/sương thiên về sáng sớm và đêm.

Trực quan hoá bộ dữ liệu

pal_many <- c("#1f77b4","#ff7f0e","#2ca02c","#d62728","#9467bd",        
               "#8c564b","#e377c2","#7f7f7f","#bcbd22","#17becf",
               "#a55194","#637939","#e6550d","#6b6ecf","#bd9e39")
Acc_raw <- Accidents

Câu lệnh: (1) Tạo một vector chứa mã màu Hex. Đây là bảng màu tuỳ chỉnh bạn sẽ tái sử dụng cho nhiều biểu đồ. (2) Sao chép dữ liệu phân tích sang đối tượng Acc_raw để giữ nguyên bản gốc.

Không gian và toạ độ

Mật độ tai nạn theo toạ độ và mức độ nghiêm trọng

ggplot(Acc_raw, aes(x = Start_Lng, y = Start_Lat)) +
    stat_bin_2d(bins = 80, aes(fill = after_stat(count))) +
    scale_fill_viridis_c(name = "Đếm") +
    facet_wrap(~Severity, nrow = 1) +
    labs(
        title = "TAI NẠN THEO TOẠ ĐỘ VÀ MỨC ĐỘ NGHIÊM TRỌNG",
        subtitle = "Phân tích mật độ qua bản đồ nhiệt 2D",
        x = "Kinh độ",
        y = "Vĩ độ" ) +
    theme_minimal() +
    theme(
        plot.title = element_text(face = "bold", hjust = 0.5),
        axis.text.x = element_text(angle = 45, hjust = 1, size = 8),
        legend.title = element_text(face = "bold"),
        panel.grid.minor = element_blank())

Câu lệnh: (1) Khởi tạo biểu đồ, lấy kinh độ làm trục X và vĩ độ làm trục Y. (2) Tính mật độ điểm theo ô lưới 2D (80 ô mỗi chiều) và dùng số điểm trong ô làm giá trị màu.(3) Áp thang màu liên tục và đặt nhãn chú giải là “Đếm”. (4) Tách biểu đồ thành nhiều ô theo mức độ nghiêm trọng, sắp trên 1 hàng. (5) Mở phần khai báo nhãn. (6-9) Đặt tên cho các trục và tiêu đề của biểu đồ. (10) Dùng nền tối giản để biểu đồ gọn và rõ. (11) Mở phần tùy biến giao diện. (12) In đậm tiêu đề và canh giữa. (13) Xoay nhãn trục X 45°, canh phải, cỡ chữ nhỏ để tránh chồng chữ. (14-15) In đậm tiêu đề chú giải và ấn lưới phụ để nền sạch hơn, tập trung vào mảng màu mật độ.

Nhận xét: Về phân bố không gian: Các vùng sáng đậm (mật độ cao ~60.000 vụ, theo thang màu) tập trung rõ ở kinh độ khoảng –120 đến –70, vĩ độ 30–45, tương ứng với các bang đông dân như California, Texas, Florida và khu vực Bờ Đông. Vùng giữa bản đồ (trung tâm nước Mỹ) mật độ nhạt hơn — phản ánh lưu lượng xe ít, khu dân cư thưa.

Về mức độ nghiêm trọng: Cả bốn cấp (1–4) có hình phân bố tương tự, cho thấy tai nạn nhẹ và nặng cùng xảy ra ở các khu vực giao thông dày đặc. Tuy nhiên, mức Severity 1–2 (tai nạn nhẹ) có mật độ dày đặc hơn hẳn màu vàng lan rộng, trong khi Severity 4 (nghiêm trọng nhất) chỉ xuất hiện lác đác tại các cụm đô thị chính.

=> Kết quả cho thấy tai nạn tập trung ở trung tâm giao thông đô thị lớn, nơi hoạt động vận tải, thương mại và du lịch mạnh, dẫn đến chi phí bảo hiểm, tắc nghẽn và thiệt hại tài sản cao. Việc này nhấn mạnh nhu cầu đầu tư hạ tầng thông minh, kiểm soát tốc độ, và giám sát hành vi lái xe tại khu vực mật độ cao (~60.000 vụ/năm) để giảm tổn thất kinh tế và nâng cao an toàn giao thông.

Mật độ tai nạn theo toạ độ

US_use <- Acc_raw %>%
  filter(!is.na(Start_Lng), !is.na(Start_Lat))
usa_states <- ggplot2::map_data("state")
bins_hex <- 80  
ggplot() +
  geom_hex(data = US_use, aes(Start_Lng, Start_Lat), 
           bins = bins_hex) +
  scale_fill_viridis_c(option = "C", trans = "log10",                             name = "Số vụ (log)") +
  geom_path(data = usa_states,                                                              aes(long, lat, group = group),                                                  color = "white", linewidth = 0.3, alpha = 0.6) +
  coord_map("albers", lat0 = 39, lat1 = 45, clip = "off") +
  labs(
    title    = "BẢN ĐỒ ĐIỂM NÓNG TAI NẠN TRÊN TOÀN NƯỚC MỸ",
    subtitle = "Hexbin thang màu log10",
    x = "Kinh độ", y = "Vĩ độ") +
  theme_void() +
  theme(
    plot.title = element_text(face = "bold", size = 12,                                                       hjust = 0.5),
    plot.subtitle = element_text(size = 12, hjust = 0.5,                                                      color = "grey30"),
    legend.position = "bottom",
    legend.key.width = grid::unit(2.6, "cm"),
    plot.margin = margin(8, 10, 8, 10))

Câu lệnh: (1–2) Lọc dữ liệu đầu vào: dùng bảng tai nạn gốc, giữ các bản ghi có kinh độ/vĩ độ hợp lệ (khác NA). (3) Nạp biên giới các bang của Mỹ từ bộ bản đồ tích hợp để vẽ đường ranh giới chồng lên. (4) Đặt tham số số “tổ” lục giác cho binning để điều khiển độ mịn của lưới. (5) Khởi tạo đối tượng đồ thị rỗng. (6–7) Thêm lớp hexbin mật độ: dùng kinh độ–vĩ độ làm tọa độ, gom theo lục giác với số bin đã đặt. Mỗi ô lục giác chứa đếm số vụ (sẽ dùng ở thang màu). (8) Thiết lập thang màu liên tục Viridis (pa-lette “C”) với chuyển đổi log10 cho giá trị đếm → làm nổi bật “điểm nóng” ở nơi có chênh lệch rất lớn. (9) Chồng lớp biên giới bang (đường ranh). (10) Đặt phép chiếu bản đồ Albers, tinh chỉnh hai vĩ tuyến chuẩn lat0, lat1. (11–14) Đặt tiêu tên. (15) Để bỏ lưới–trục mặc định, tập trung vào lớp không gian. (16–21) Canh chỉnh theme.

Nhận xét: Các vùng màu vàng–cam sáng (mật độ cao) tập trung dày đặc ở bờ Đông, đặc biệt quanh New York, Washington D.C., Florida và Đông Nam Texas – California. Các khu vực Trung Tây và miền Tây Bắc có mật độ thưa (màu tím đậm), cho thấy lưu lượng xe và dân cư thấp hơn. => Bản đồ thể hiện mật độ tai nạn theo không gian, giúp nhận diện vùng rủi ro. Việc dùng log10 giúp tránh “cháy màu” ở đô thị nhưng vẫn làm rõ các tuyến đường thưa dân; có thể tinh chỉnh bins (60–100) hoặc độ trong suốt (alpha) để tối ưu hiển thị.

Toạ độ tai nạn theo ca và tầm nhìn

df <- Acc_raw %>%
  mutate( vis_raw = suppressWarnings(as.numeric(`Visibility(mi)`))) %>%
  filter(!is.na(vis_raw)) %>%
  mutate(
    q1  = quantile(vis_raw, 0.01, na.rm = TRUE),
    q99 = quantile(vis_raw, 0.99, na.rm = TRUE),
    vis = pmin(pmax(vis_raw, q1), q99)   )
ggplot() +
  geom_polygon(data = ggplot2::map_data("state"),
               aes(long, lat, group = group),
               fill = "#f3f4f6", color = "white", linewidth = 0.25) +
  geom_point(data = df,
             aes(Start_Lng, Start_Lat, fill = time_period, size = vis),
             shape = 21, color = "white", stroke = 0.1, alpha = 0.65) +
  facet_wrap(~ time_period, ncol = 2) +
  scale_size_continuous(
    range  = c(0.6, 3.4),         
    trans  = "sqrt",               
    breaks = c(0, 10, 25, 50, 75, 100),
    limits = c(0, 100),          
    name   = "Tầm nhìn (mi)") +
  scale_fill_manual(values = c("Đêm"="#3B5BFE",
        "Sáng"="#FF8C00","Chiều"="#2ECC71","Tối"="#E74C3C"),
                    name = "Ca") +
  coord_map("albers", lat0 = 39, lat1 = 45, clip = "on") +
  guides(
    fill = guide_legend(override.aes = list(size = 3,
                                            alpha = 0.9, stroke = 0.2)),
    size = guide_legend(order = 1)) +
  labs(
    title = "TỌA ĐỘ TAI NẠN THEO CA & TẦM NHÌN",
    subtitle = "Kích thước điểm tầm nhìn (mi)",   
    x = "Kinh độ", y = "Vĩ độ") +
  theme_minimal(base_size = 11) +
  theme(
    panel.grid.minor = element_blank(),
    panel.grid.major = element_line(color = "#e6e8ef", linewidth = 0.25),
    strip.background = element_rect(fill = "#eef0f5", color = NA),
    strip.text = element_text(face = "bold"))

Câu lệnh: (1–3) Tạo biến vis_raw: chuyển cột Visibility(mi) sang dạng số, loại bỏ giá trị thiếu (NA). (4–5) Tính ngưỡng ngoại lệ: q1 và q99 lần lượt là phân vị 1% và 99% của tầm nhìn, giúp cắt bỏ các giá trị cực đoan. (6) Giới hạn tầm nhìn hợp lý: ép tất cả giá trị nằm trong khoảng [q1, q99] để giảm nhiễu. (7–10) Khởi tạo ggplot bản đồ. (11–13) Thêm điểm tai nạn: mỗi điểm biểu thị vị trí bắt đầu (Start_Lng, Start_Lat), tô màu theo thời điểm (time_period), kích thước theo tầm nhìn (vis), viền trắng và độ trong suốt 0.65. (14) Chia ô theo thời điểm. (15–20) Chuẩn hóa kích thước điểm: ánh xạ kích thước theo vis từ 0.6–3.4, thang căn bậc hai (sqrt), giới hạn 0–100 mi, nhãn “Tầm nhìn (mi)”. (21–23) Đặt bảng màu thời điểm trong ngày. (24) Căn bản đồ kiểu Albers: điều chỉnh phép chiếu bản đồ với vĩ tuyến 39–45, bật clip = “on” để không tràn biên. (25–26) Tinh chỉnh chú giải. (27–32) Thêm tiêu đề, nhãn trục. (33–34) Đặt giao diện tối giản. (35–39) Tùy chỉnh chi tiết giao diện.

Nhận xét:

Biểu đồ thể hiện phân bố tai nạn theo ca và tầm nhìn (mi) trên bản đồ Hoa Kỳ. Kích thước điểm càng lớn → tầm nhìn càng xa.

Ca Sáng (màu cam) và Chiều (màu xanh lá) có mật độ tai nạn dày đặc, phủ gần như toàn quốc, phản ánh lưu lượng phương tiện cao vào hai khung giờ cao điểm trong ngày.

Ca Tối (màu đỏ) vẫn ghi nhận nhiều vụ, song kích thước điểm nhỏ hơn, cho thấy tầm nhìn bị hạn chế hơn — yếu tố rủi ro điển hình khi ánh sáng yếu.

Ca Đêm (màu xanh dương) ít vụ hơn rõ rệt, tập trung dọc theo các tuyến cao tốc hoặc khu đô thị lớn, tầm nhìn thấp, thể hiện ảnh hưởng của bóng tối, mệt mỏi và điều kiện chiếu sáng yếu.

Nhìn chung, tai nạn tập trung cao ở Sáng và Chiều, khi hoạt động giao thông sôi động nhất; trong khi Tối và Đêm có ít vụ hơn nhưng rủi ro cao hơn do tầm nhìn giảm. Điều này gợi ý cần tăng cường chiếu sáng, biển cảnh báo và tuần tra ban đêm, đồng thời điều tiết giao thông hợp lý trong giờ cao điểm ban ngày.

Toạ độ theo ca trong ngày

ggplot(Acc_raw, aes(Start_Lng, Start_Lat, color=time_period)) +
  geom_point(alpha=.08, size=.5) +
  stat_density_2d(aes(linetype=time_period), color="black", bins=5)+
  scale_color_manual(values=pal_many) +
  labs(title="TOẠ ĐỘ THEO CA TRONG NGÀY (CONTOURS)") +
  theme_minimal()

Câu lệnh: (1) Khởi tạo biểu đồ, trục x và y lần lượt là kinh độ và vĩ độ, màu thể hiện ca trong ngày. (2) Thêm các điểm tai nạn, hiển thị mờ và nhỏ để giảm chồng lấn. (3) Vẽ các đường mật độ thể hiện khu vực tập trung điểm, đường viền đen, phân tách thành 5 mức, kiểu đường thay đổi theo ca. (4) Áp dụng bảng màu thủ công cho từng ca trong ngày. (5) Thêm tiêu đề cho biểu đồ, mô tả tọa độ tai nạn theo ca trong ngày. (6) Dùng giao diện tối giản để làm nổi bật vùng mật độ.

Nhận xét:

Các vùng đậm nét tập trung rõ ở California (Los Angeles–San Francisco), Florida (Miami–Orlando), Texas (Dallas–Houston) và Bờ Đông (New York–Washington D.C.) — đều là khu vực đô thị lớn, lưu lượng phương tiện cao.

Đường viền của ca Chiều (nét gạch ngắn) và Sáng (nét chấm) chiếm ưu thế, cho thấy tai nạn tập trung vào hai khung giờ cao điểm, trùng với kết quả từ biểu đồ thời gian.

Ca Tối và Đêm có ít vùng viền hơn, diện tích nhỏ, phản ánh lượng tai nạn thấp hơn nhưng vẫn tập trung ở đô thị chính, đặc biệt quanh các tuyến cao tốc và cảng.

Nhìn chung, biểu đồ cho thấy tai nạn giao thông có xu hướng tập trung tại các trung tâm đô thị lớn và hành lang giao thông liên vùng, với đỉnh rủi ro vào ca Sáng và Chiều.

Bản đồ điểm nóng theo bang (Top 4 State)

top4 <- names(sort(table(Acc_raw$State), TRUE))[1:4]
map_us <- map_data("state") |>
  mutate(State = state.abb[match(region, tolower(state.name))]) |>
  filter(State %in% top4)

Acc_raw |> filter(State %in% top4) |>
  ggplot(aes(Start_Lng, Start_Lat)) +
  stat_bin_hex(bins = 55, alpha = 0.9) +
  geom_path(data = map_us, aes(long, lat, group = group),
            color = "white", linewidth = 0.3) +
  facet_wrap(~State, ncol = 4, scales = "free") +
  scale_fill_viridis_c(option = "B", trans = "log10",
                       name = "Số vụ (log10)") +
  coord_quickmap() +
  labs(
    title = "BẢN ĐỒ ĐIỂM NÓNG THEO BANG (TOP 4)",
    subtitle = "Hexbin, thang log10, lưới kinh–vĩ độ",
    x = "Kinh độ", y = "Vĩ độ"
  ) +
  theme_minimal(base_size = 11) +
  theme(
    axis.text.x = element_text(angle = 35, hjust = 1, size = 8),
    strip.text = element_text(face = "bold"),
    plot.title = element_text(face = "bold"))

Câu lệnh: (1–2) Lấy 4 bang có nhiều tai nạn nhất. (3–4) Tạo bản đồ ranh giới các bang của Mỹ. (5–6) Lọc dữ liệu gốc để chỉ còn các bản ghi trong 4 bang đó. (7–8) Khởi tạo biểu đồ. (9–10) Vẽ đường viền ranh giới bang màu trắng, nét mảnh để tách biệt các bang. (11) Chia biểu đồ thành 4 ô nhỏ, mỗi ô là một bang. (12–13) Dùng thang màu log10 giúp thể hiện rõ chênh lệch giữa vùng ít và nhiều tai nạn, thêm chú giải. (14) Dùng bản đồ nhanh để tỷ lệ giữa trục x–y chuẩn hơn với bản đồ thực. (15–18) Thêm tiêu đề, phụ đề và nhãn trục. (19–20) Áp dụng theme tối giản để làm nổi bật các ô hexbin. (21–24) Tùy chỉnh giao diện: x nghiên 35°.

Nhận xét:

Bản đồ điểm nóng theo bang (Top 4) cho thấy sự tập trung rõ rệt của các vụ tai nạn giao thông tại những khu vực đô thị đông dân và có lưu lượng phương tiện lớn.

Ở California (CA), mật độ tai nạn cao nhất xuất hiện dọc ven biển, đặc biệt quanh các thành phố lớn như Los Angeles, San Francisco và San Diego – nơi hạ tầng giao thông dày đặc và mật độ xe cộ cao.

Florida (FL) có các điểm nóng trải dài dọc bờ biển phía Đông, tập trung ở vùng Nam Florida (Miami, Fort Lauderdale) và khu vực trung tâm (Orlando, Tampa). Điều này phản ánh ảnh hưởng của mạng lưới đường ven biển, du lịch phát triển và dân cư đông đúc.

Tại Texas (TX), các cụm mật độ cao tập trung tại Houston, Dallas–Fort Worth, Austin và San Antonio. Phân bố có phần rải rác hơn do diện tích bang lớn, nhưng vẫn xoay quanh các đô thị trọng điểm.

Trong khi đó, Virginia (VA) có mật độ tai nạn cao nhất ở khu vực Đông Bắc gần Washington D.C., giảm dần về phía Tây và Nam – thể hiện sự khác biệt giữa vùng đô thị và nông thôn.

Tổng thể, bản đồ cho thấy tai nạn giao thông ở Mỹ mang tính đô thị hóa mạnh, tập trung chủ yếu tại các trung tâm kinh tế – giao thông lớn. Việc sử dụng thang log10 giúp làm nổi bật sự chênh lệch mật độ giữa các khu vực đông và thưa tai nạn mà vẫn giữ được độ tương phản rõ ràng, tránh hiện tượng “cháy màu” trên biểu đồ.

Thời gian

Phân bố số vụ tai nạn theo giờ và ca trong ngày

peaks <- Acc_raw |> filter(!is.na(Hour), !is.na(time_period)) |>
  count(time_period, Hour, name = "n") |> group_by(time_period) |>
  slice_max(order_by = n, n = 1) |> ungroup()
base <- Acc_raw |> filter(!is.na(Hour), !is.na(time_period))
ggplot(base, aes(x = Hour, y = time_period, fill = after_stat(x))) +
  ggridges::geom_density_ridges_gradient(scale = 1.5, rel_min_height = 0.01,
                                         color = "grey20", size = .25, bandwidth = .55) +
  geom_vline(data = peaks, aes(xintercept = Hour), color = "#505050",
             linetype = "22", linewidth = .35) +
  geom_label(data = peaks, aes(x = Hour, y = time_period,
             label = paste0("Peak: ", Hour, "h")),
             fill = "white", color = "#222", size = 3.2, vjust = -0.45) +
  scale_x_continuous(breaks = seq(0, 23, 2), limits = c(0, 23)) +
  scale_fill_viridis_c(option = "M", direction = -1, guide = "none") +
  labs(title = "PHÂN BỐ GIỜ TAI NẠN THEO CA", x = "Giờ trong ngày", y = NULL) +
  theme_minimal(base_size = 13) +
  theme(plot.title = element_text(face="bold", size=14),
        axis.text.y = element_text(face="bold"))

Câu lệnh: (1–3) Đếm số vụ theo từng ca và giờ, rồi chọn mỗi ca 1 giờ có số vụ lớn nhất (giờ “peak”). (4) Tạo dữ liệu vẽ: chỉ giữ các bản ghi có giờ và ca hợp lệ. (5) Khởi tạo biểu đồ: trục ngang là giờ, trục dọc là ca; màu tô của dải sẽ đổi theo giá trị giờ. (6–7) Vẽ các dải ridge có tô gradient theo giờ; tăng khoảng cách giữa các dải, bỏ đuôi rất mỏng, viền xám mảnh; làm mượt bằng băng thông 0,55. (8–9) Kẻ vạch dọc tại giờ peak của từng ca, màu xám, nét đứt mảnh để dễ đối chiếu. (10–12) Gắn nhãn “Peak”. (13) Đặt vạch chia trục X mỗi 2 giờ trong khoảng 0–23. (14) Dùng thang màu Viridis biến thiên theo giờ, ẩn chú giải để gọn. (15) Thêm tiêu đề. (16) Áp dụng kiểu giao diện tối giản.m(17–18) In đậm tiêu đề và nhãn trục.

Nhận xét: Ca sáng đạt đỉnh quanh 7h, mật độ cao và tập trung, trùng với giờ đi làm, đi học. Ca chiều đỉnh 16h, kéo dài từ 14–18h, phản ánh cao điểm tan tầm. Ca tối đỉnh 19h, lan đến 22h, gắn với hoạt động giải trí và tầm nhìn kém. Ca đêm tuy ít nhưng vẫn có đỉnh nhỏ khoảng 6h sáng, do mệt mỏi cuối ca.

=> Tai nạn tập trung vào ba khung 6–8h, 16–18h và 19–21h, hàm ý cần tăng cường giám sát, chiếu sáng và điều tiết giao thông trong các giai đoạn này.

Xu hướng theo ngày và ca

daily <- Acc_raw %>% count(Date, time_period, name = "n")
ggplot(daily, aes(Date, n, color=time_period, group=time_period)) +
  geom_line(alpha=.4) +
  geom_smooth(se=FALSE, method="gam", formula = y ~ s(x, k = 20)) +
  scale_color_manual(values=pal_many) +
  labs(title="SỐ VỤ TAI NẠN THEO NGÀY VÀ CA", x="Ngày", y="Số vụ") +
  theme_minimal()

Câu lệnh: (1) Gom số liệu theo ngày (Date) và ca trong ngày (time_period), đếm số vụ vào cột n. (2) Tạo khung vẽ, trục x là ngày, trục y là số vụ, tô màu theo ca. (3) Vẽ đường thời gian hằng ngày cho từng ca. (4) Thêm đường xu thế trơn riêng cho từng ca dùng GAM, k=20 nút giúp nắm xu hướng/chu kỳ mà không bám sát nhiễu ngày-ngày. (5) Dùng bảng màu tùy chỉnh cho các ca. (6-7) Đặt tiêu đề nhãn trục và dùng theme tối giản.

Nhận xét: Biểu đồ cho thấy ca Chiều (xanh lá) luôn chiếm tỷ trọng tai nạn cao nhất, kế đến là Sáng, còn Tối và Đêm thấp hơn đáng kể → trùng với quy luật giờ cao điểm sáng và chiều.

Về tính mùa vụ, cả Sáng và Chiều tăng mạnh vào Q1–Q2 (tháng 3–5), giảm vào hè (Q3) rồi nhích lại Q4, phù hợp với nhịp học tập, lao động và thương mại cuối năm, đồng thời chịu ảnh hưởng thời tiết chuyển mùa như mưa xuân, sương mù. Biến động ngày–ngày dao động mạnh quanh xu hướng, nhất là ở ca Chiều, phản ánh ảnh hưởng của cuối tuần, lễ, hoặc thời tiết cực đoan.

Phân bố số vụ theo giờ và ngày trong tuần

bubble <- Acc_raw |> 
  filter(!is.na(Hour), !is.na(Day)) |>
  mutate(Day = factor(Day, levels = c("Monday","Tuesday","Wednesday","Thursday","Friday","Saturday","Sunday"))) |>
  count(Day, Hour, name = "n")
ggplot(bubble, aes(Hour, Day)) +
  geom_point(aes(size = n, fill = n), shape = 21, colour = alpha("grey10", .25), stroke = .35) +
  scale_fill_scico(palette = "lajolla", labels = label_number(big.mark=".", decimal.mark=","), name = "Số vụ (màu)") +
  scale_size_area(max_size = 10, breaks = c(10000, 20000),
                  labels = label_number(big.mark=".", decimal.mark=","), name = "Số vụ") +
  scale_x_continuous(breaks = seq(0, 23, 2), minor_breaks = 0:23) +
  guides(size = guide_legend(override.aes = list(shape = 21, fill = "grey60", colour = NA))) +
  labs(title = "PHÂN BỐ SỐ VỤ NGÀY THEO GIỜ",
       subtitle = "Màu và kích thước cùng biểu diễn số vụ",
       x = "Giờ", y = "Ngày trong tuần") +
  theme_minimal(base_size = 10) +
  theme(panel.grid.major.y = element_blank(),
        panel.grid.minor = element_blank(),
        panel.grid.major.x = element_line(colour = alpha("grey35", .14), linewidth = .7),
        legend.position = "right",
        legend.title = element_text(face = "bold"),
        plot.title = element_text(face = "bold", size = 10),
        plot.subtitle = element_text(colour = "grey40"))

Câu lệnh:(1–3) Lọc dữ liệu chỉ giữ bản ghi có giờ và ngày hợp lệ, đồng thời sắp xếp thứ tự các ngày trong tuần từ Thứ Hai → Chủ Nhật. (4) Đếm số vụ theo từng giờ trong từng ngày → chuẩn bị cho biểu đồ bong bóng. (5–7) Vẽ biểu đồ: trục x = giờ, trục y = ngày trong tuần, kích thước và màu bóng thể hiện số vụ. (8–10) Thiết lập kích thước bóng tối đa và thang chia giá trị, định dạng số kiểu Việt. (11) Chia mốc trục x cách 2 giờ/lần để dễ đọc. (12) Tùy chỉnh chú giải: mẫu minh họa hình tròn tô xám, đặt bên phải. (13–15) Tiêu đề, phụ đề và nhãn trục: giải thích ý nghĩa màu và kích thước bóng.(16–22) Dùng theme minimal: ẩn lưới ngang, giữ lưới dọc mờ để canh giờ, căn chỉnh tiêu đề/chú giải đậm rõ, phông chữ gọn gàng.

Nhận xét: Tai nạn tập trung rõ vào giờ cao điểm 7–9h sáng và 15–18h chiều, lặp lại từ thứ Hai đến thứ Sáu — phản ánh nhịp đi lại giờ đi làm và tan ca. Thứ Sáu có số vụ cao nhất, với bong bóng sáng và to ở cả hai đỉnh sáng–chiều → giai đoạn giao thông dày đặc trước cuối tuần. Cuối tuần (Saturday–Sunday) tai nạn giảm đáng kể, bong bóng nhỏ và sẫm hơn, do lưu lượng đi lại giảm và ít giao thông giờ cao điểm.

Phân bố số vụ theo 12 tháng

d <- Acc_raw |> count(Month = month(Date, label=TRUE, abbr=TRUE),  Hour, name="n")
avg <- d |> group_by(Hour) |> summarise(avg=mean(n), .groups="drop")
pk  <- d |> group_by(Month) |> slice_max(n, n=1) |> ungroup()
Ymax <- max(d$n); bands <- data.frame(xmin=c(7,15), xmax=c(9,18), ymin=-Inf, ymax=Inf)
pk <- pk |> mutate(lbl=paste0(Hour,"h\n",comma(n,big.mark=".",decimal.mark=",")),
                   y=n+0.08*Ymax)
ggplot(d,aes(Hour,n))+
  geom_line(data=avg,aes(Hour,avg),inherit.aes=F,col="#9AA0A6",linetype=22, linewidth=1)+
  geom_line(col="#1F4B99",linewidth=1.2)+
  geom_point(data=pk,col="#C43E6D",size=2.8)+
  geom_text_repel(data=pk,aes(y=y,label=lbl),col="#C43E6D",size=3.2,
          nudge_y=0.03*Ymax,segment.colour="#C43E6D",min.segment.length=0,seed=1)+
  facet_wrap(~Month,ncol=4)+                            
  scale_x_continuous(breaks=c(0,6,12,18,23))+
  scale_y_continuous(labels=label_number(big.mark=".",decimal.mark=","),
                     expand=expansion(mult=c(0.02,0.08)))+
  labs(title="SỐ VỤ THEO GIỜ – 12 THÁNG",
       subtitle="· Đường xám: trung bình năm · Dấu đỏ: giờ đỉnh",
       x="Giờ",y="Số vụ")+
  coord_cartesian(clip="off")+
  theme_minimal(base_family="Times New Roman",base_size=14)+
  theme(panel.grid.minor=element_blank(),
        strip.text=element_text(face="bold", size=12),
        axis.text=element_text(size=11),
        axis.title=element_text(size=13),
        plot.title=element_text(face="bold",hjust=.5, size=20),
        plot.subtitle=element_text(hjust=.5, size=12),
        panel.spacing=unit(10,"pt"),
        plot.margin=margin(10,20,10,10))

Câu lệnh: (1–3) Tạo bảng đếm số vụ theo tháng và giờ, tính trung bình theo giờ và lấy giờ đỉnh mỗi tháng. (4–5) Xác định giá trị cao nhất để đặt nhãn, rồi tạo nhãn hiển thị giờ và số vụ, dời nhẹ lên trên để tránh chồng. (6–10) Vẽ biểu đồ: đường xám là trung bình năm, đường xanh là theo tháng, dấu hồng là giờ đỉnh, và nhãn hồng hiển thị giá trị tại điểm cao nhất. (13–15) Chia ô theo 12 tháng, trục x gồm các mốc chính (0–23h), trục y định dạng kiểu Việt và chừa khoảng trống cho nhãn. (16–18) Thêm tiêu đề, phụ đề giải thích ký hiệu, nhãn trục “Giờ” và “Số vụ”. (19–28) Dùng theme minimal (Times New Roman), ẩn lưới phụ, chữ đậm dễ đọc, căn giữa tiêu đề, nới khoảng cách và lề cho bố cục gọn gàng.

Nhận xét: Biểu đồ “Số vụ theo giờ – 12 tháng” cho thấy quy luật rõ rệt về thời gian xảy ra tai nạn trong ngày và mức độ biến động theo tháng.

Cụ thể, đỉnh tai nạn thường xuất hiện vào 15h–17h, trùng với khung giờ tan tầm chiều — khi lưu lượng giao thông lớn và dễ xảy ra va chạm. Ngoài ra, nhiều tháng cũng xuất hiện đỉnh phụ nhỏ vào khoảng 7h–9h sáng, là giờ cao điểm buổi sáng.

Nhìn chung, các tháng có hình dạng đường cong khá tương đồng, chỉ khác nhau về cường độ. Các tháng mùa khô (từ tháng 4 đến 8) có xu hướng cao hơn trung bình năm, có thể do thời tiết thuận lợi khiến tần suất di chuyển tăng. Đường xám thể hiện mức trung bình năm, giúp so sánh sự khác biệt giữa từng tháng; còn dấu đỏ làm nổi bật giờ cao điểm nhất của mỗi tháng, cho phép nhận diện nhanh điểm rủi ro về thời gian. Tổng thể, biểu đồ minh họa chu kỳ rõ rệt theo giờ trong ngày, phản ánh hành vi di chuyển mang tính thói quen và gợi ý thời điểm nên tăng cường kiểm soát giao thông để giảm thiểu rủi ro.

Phân bố Severity theo giờ và ca

df_hour <- Acc_raw |>
  filter(!is.na(Hour), !is.na(time_period)) |>
  mutate(time_period = factor(time_period, levels = c("Đêm","Sáng","Chiều","Tối"))) |>
  group_by(time_period, Hour) |>
  summarise(sev = mean(Severity, na.rm = TRUE), .groups = "drop") |>
  group_by(time_period) |>
  mutate(is_max = sev == max(sev, na.rm = TRUE)) |>
  ungroup()
pal <- c("Đêm"="#6DAA2C","Sáng"="#23B5D3","Chiều"="#F28D85","Tối"="#8B7CF6") 
ggplot(df_hour, aes(Hour, sev, group = time_period, color = time_period)) +
  geom_smooth(aes(fill = time_period), method = "gam", formula = y ~ s(x, k = 6),
              se = TRUE, alpha = .12, color = NA) +
  geom_line(linewidth = 1.1) +
  geom_point(shape = 21, stroke = 1.1, fill = "white", size = 2.6) +
  ggrepel::geom_label_repel(
    data = subset(df_hour, is_max),
    aes(label = paste0("Max: ", scales::number(sev, accuracy = 0.001, decimal.mark = ","))),
    seed = 4, size = 3.2, label.size = 0, fill = scales::alpha("white", .85),
    min.segment.length = .1) +
  scale_color_manual(values = pal, name = "Ca") +
  scale_fill_manual(values = pal, guide = "none") +
  scale_x_continuous(breaks = seq(0, 23, 2)) +
  scale_y_continuous(labels = label_number(accuracy = 0.001, decimal.mark = ",")) +
  labs(title = "MỨC ĐỘ NGHIÊM TRỌNG THEO GIỜ VÀ CA",
       subtitle = "Đường màu theo ca; vùng mờ là xu hướng GAM (±SE)",
       x = "Giờ", y = "Severity (trung bình)") +
  theme_minimal(base_size = 13) +
  theme(plot.title = element_text(face = "bold", size = 16),
        legend.title = element_text(face = "bold"),
        panel.grid.minor = element_blank())

Câu lệnh: (1–4) Lọc dữ liệu có giờ và ca hợp lệ, sắp xếp theo thứ tự Đêm–Sáng–Chiều–Tối, rồi tính mức độ nghiêm trọng trung bình (Severity) cho từng giờ trong từng ca. (5) Gắn cờ is_max = TRUE ở giờ có mức nghiêm trọng cao nhất trong mỗi ca để gắn nhãn sau này. (6) Tạo bảng màu cho 4 ca. (7–8) Khởi tạo biểu đồ: trục hoành là giờ, trục tung là mức độ nghiêm trọng trung bình, màu theo ca. (9) Thêm đường xu hướng GAM có vùng mờ ±SE để biểu diễn biến động mượt. (10–11) Vẽ đường chính và điểm trắng thể hiện giá trị trung bình từng giờ. (12–17) Thêm nhãn “Max” cho giờ có mức nghiêm trọng cao nhất mỗi ca. (18–20) Tùy chỉnh thang màu, trục x (0–23h), trục y hiển thị theo định dạng số Việt. (21–26) Đặt tiêu đề, phụ đề, nhãn trục. (27–30) Giao diện tối giản: tiêu đề và chú giải in đậm, bỏ lưới phụ.

Nhận xét: Ca Đêm (xanh lá) có mức nghiêm trọng cao nhất, đỉnh khoảng 2,117 vào 1–2h sáng, sau đó giảm dần khi gần sáng. Ca Tối (tím) đứng thứ hai, đạt đỉnh 2,118 quanh 20–21h, cho thấy rủi ro cao vào ban đêm khi tầm nhìn kém. Ca Sáng (xanh dương) và Chiều (cam) có mức thấp hơn (≈2,08–2,09), ổn định hơn nhưng vẫn có đỉnh nhẹ vào giờ cao điểm.

=> Nhìn chung, tai nạn ban đêm và buổi tối ít hơn về số lượng nhưng nghiêm trọng hơn, do yếu tố mệt mỏi, tốc độ cao và ánh sáng yếu. Hàm ý: cần tăng chiếu sáng, kiểm soát tốc độ và tuần tra đêm để giảm rủi ro nặng.

Mật độ Severity theo ca trong ngày

Acc_raw %>%
  filter(!is.na(Severity), !is.na(time_period)) %>%
  mutate(time_period = factor(time_period,
         levels = c("Đêm","Sáng","Chiều","Tối"))) %>%
  ggplot(aes(x = Severity, y = time_period)) +
  ggdist::stat_slab(aes(thickness = after_stat(pdf), fill = after_stat(pdf)),
                    normalize = "groups", adjust = 0.6,
                    colour = "grey20", linewidth = 0.25) +
  ggdist::stat_pointinterval(aes(x = Severity),
                             point_interval = ggdist::mean_qi,
                             .width = c(0.5, 0.8, 0.95),
                             position = position_nudge(y = 0.07),
                             size = 2, colour = "grey15") +
  scale_fill_scico(palette = "lajolla", direction = 1, name = "Mật độ") +
  labs(title = "Mật độ Severity theo ca trong ngày",
       subtitle = "Độ dày & màu biểu diễn mật độ; điểm/khoảng = mean & QI",
       x = "Severity", y = "Ca") +
  theme_minimal(base_size = 13) +
  theme(plot.title = element_text(face = "bold"),
        panel.grid.minor = element_blank())

Câu lệnh: (1–3) Lọc dữ liệu có Severity và ca (time_period) hợp lệ, sắp xếp theo thứ tự “Đêm–Sáng–Chiều–Tối”. (4–5) Tạo biểu đồ với trục x = Severity, trục y = Ca. (6–8) Vẽ dải mật độ (slab): độ dày và màu thể hiện mật độ xác suất; mỗi ca chuẩn hoá riêng. (9–12) Thêm trung bình và khoảng tin cậy (50%, 80%, 95%) bằng chấm và thanh ngang. (13) Thang màu scico thể hiện mật độ. (14–17) Gắn tiêu đề, phụ đề và nhãn trục. (18–20) Dùng theme tối giản, in đậm tiêu đề, bỏ lưới phụ.

Nhận xét: Các dải mật độ (slab) rất hẹp quanh mức Severity ≈ 2, nghĩa là phần lớn các vụ tai nạn có mức nghiêm trọng trung bình thấp, ít dao động. Độ dày và màu của dải gần như tương đương giữa các ca → không có sự khác biệt lớn giữa các khung giờ về độ nghiêm trọng trung bình. Một vài điểm cực trị nhỏ ở mức Severity = 3–4 cho thấy một số ít vụ nghiêm trọng cao, nhưng tần suất cực thấp (mật độ gần như bằng 0).

Phân bố nhiệt độ theo ca trong ngày

df <- Acc_raw %>%
  filter(!is.na(temp_c), !is.na(time_period)) %>%
  mutate(time_period = factor(time_period, levels = c("Đêm","Sáng","Chiều","Tối")))
ggplot(df, aes(x = time_period, y = temp_c, fill = time_period)) +
  gghalves::geom_half_violin(side = "l", width = .95, trim = TRUE,
                             alpha = .85, colour = NA) +
  geom_boxplot(width = .14, outlier.alpha = 0, fill = "white",
               colour = "grey25", linewidth = .35) +
  ggdist::stat_dots(side = "right",  
                    scale = 0.6,     
                    dotsize = 1.1,
                    alpha = .25, colour = "grey30") +
  scale_fill_scico_d(palette = "lajolla", direction = 1, guide = "none") +
  coord_flip() +
  labs(title = "PHÂN BỐ NHIỆT ĐỘ THEO CA - RAINCLOUD",
       x = NULL, y = "Nhiệt độ (°C)") +
  theme_minimal(base_size = 13) +
  theme(panel.grid.minor = element_blank())

Câu lệnh: (1–3) Lọc các bản ghi có nhiệt độ & ca hợp lệ, sắp xếp ca theo thứ tự Đêm → Sáng → Chiều → Tối. (4) Thiết lập trục: ngang = Ca, dọc = Nhiệt độ (°C), tô màu theo ca. (5) Vẽ nửa violin bên trái để thấy phân bố nhiệt độ của từng ca. (7–8) Thêm boxplot mảnh (ẩn outlier) ở giữa để đọc trung vị & tứ phân vị rõ ràng. (9–12) Rải chấm dữ liệu bên phải (dotplot) để hiển thị điểm quan sát thô. (13) Dùng bảng màu đồng nhất theo ca, tắt chú giải cho gọn. (14) Xoay ngang toàn bộ để nhãn ca dễ đọc, so sánh nhanh. (15–17) Gắn tiêu đề, nhãn trục (ẩn trục x vì đã là tên ca). (18–19) Áp theme tối giản, bỏ lưới phụ để bố cục sạch.

Nhận xét

Đêm: Nhiệt độ phân bố chủ yếu ở mức thấp (âm hoặc gần 0°C), dao động hẹp → không khí lạnh và ổn định hơn, phản ánh đặc trưng ban đêm ít chịu tác động của bức xạ mặt trời.

Sáng: Nhiệt độ trung bình tăng nhẹ, độ phân tán lớn hơn do giai đoạn chuyển tiếp từ lạnh sang ấm; xuất hiện biến thiên mạnh giữa các khu vực.

Chiều: Là khoảng có nhiệt độ cao nhất trong ngày, trung vị dịch sang phải rõ rệt → tích tụ nhiệt tối đa sau khi mặt trời chiếu mạnh và lâu.

Tối: Nhiệt độ bắt đầu giảm dần, phân bố hẹp hơn nhưng vẫn cao hơn sáng sớm → giai đoạn giải nhiệt chậm của không khí.

Thời tiết và điều kiện môi trường

Nhiệt độ và tầm nhìn theo trạng thái thời tiết

pal_weather <- c(
  "Trời quang" = "#2A9D8F", "Nhiều mây" = "#E9C46A",
  "Mưa" = "#264653", "Tuyết" = "#E76F51",
  "Sương mù" = "#F4A261", "Khác" = "grey50")
Acc_raw %>%
  filter(!is.na(weather_simple), !is.na(`Temperature(F)`), !is.na(`Visibility(mi)`)) %>%
  mutate(weather_simple = fct_lump_n(weather_simple, 5, other_level = "Khác")) %>%
  ggplot(aes(`Temperature(F)`, `Visibility(mi)`, color = weather_simple)) +
  geom_point(alpha = 0.05, size = 1) +
  geom_smooth(method = "gam", formula = y ~ s(x, k = 5), se = TRUE, linewidth = 1.2) +
  facet_wrap(~weather_simple, ncol = 3, scales = "free") +
  scale_color_manual(values = pal_weather) +
  labs(
    title = "MỐI QUAN HỆ GIỮA NHIỆT ĐỘ VÀ TẦM NHÌN THEO TRẠNG THÁI THỜI TIẾT",
    subtitle = "Đường cong thể hiện xu hướng trung bình trong từng điều kiện thời tiết",
    x = "Nhiệt Độ (F)", y = "Tầm Nhìn (miles)") +
  theme_minimal(base_size = 10) +
  theme(
    plot.title = element_text(face = "bold", size = 10),
    plot.subtitle = element_text(color = "grey50"),
    strip.text = element_text(face = "bold", size = 10),
    panel.grid.minor = element_blank(),
    legend.position = "none")

Câu lệnh: (1–3) Tạo bảng màu cho 6 trạng thái thời tiết; lọc dữ liệu hợp lệ (Temperature, Visibility, weather_simple) và gộp nhóm hiếm vào “Khác”. (4–6) Thiết lập biểu đồ với trục x = Nhiệt độ (°F), y = Tầm nhìn (mi); màu biểu diễn trạng thái thời tiết. (7–9) Vẽ điểm mờ thể hiện phân bố thực tế; geom_smooth(method=“gam”) vẽ đường xu hướng trung bình và vùng sai số (±SE). (10–12) facet_wrap() tách từng điều kiện thời tiết thành ô riêng (3 cột), trục độc lập để dễ so sánh. (13–15) Dùng bảng màu thủ công giúp giữ màu nhất quán. (16–18) Thêm tiêu đề, phụ đề, nhãn trục. (19–23) Làm bố cục gọn, in đậm tiêu đề, ẩn lưới phụ và chú giải.

Nhận xét

Mưa: Tầm nhìn thấp, dao động chủ yếu 0–10 miles. Đường xu hướng GAM tăng mạnh khi nhiệt độ cao hơn ⇒ mưa ấm, nhẹ giúp tầm nhìn cải thiện; khi lạnh (mưa lạnh/ẩm) tầm nhìn giảm.

Trời quang: Tầm nhìn cao nhất (~100 miles), đường gần phẳng → điều kiện trong, nhiệt độ hầu như không ảnh hưởng.

Nhiều mây: Tầm nhìn trung bình (~80 miles), đường GAM hơi cong dạng vòm – tốt nhất ở nhiệt độ trung bình, giảm khi quá lạnh hoặc quá nóng.

Sương mù: Tầm nhìn rất thấp (< 5 miles), đường GAM tăng rõ theo nhiệt độ ⇒ không khí ấm giúp sương tan.

Tuyết: Tầm nhìn 0–10 miles, tốt nhất quanh 25–35°F rồi giảm khi quá lạnh hoặc quá ấm (mưa tuyết, băng giá).

Khác: Xu hướng tăng tuyến tính nhẹ, nhiệt độ cao giúp nhìn xa hơn

=> Tầm nhìn bị ảnh hưởng mạnh nhất bởi sương mù, tuyết và mưa, trong khi trời quang và nhiều mây duy trì điều kiện quan sát tốt.

Tầm nhìn theo trạng thái thời tiết

pal <- c("#4A7C59","#A63D40","#D3A556","#738678","#3F4B3F",
         "#6D7B8C","#9A9A9A","#4C72B0","#F4A261","#2A9D8F")
Acc_raw %>%
  filter(!is.na(weather_simple), !is.na(`Visibility(mi)`)) %>%
  mutate(weather_simple = fct_lump_n(weather_simple, 10, other_level = "Khác")) %>%
  filter(weather_simple != "Khác") %>%                                   
  mutate(weather_simple = reorder(weather_simple, `Visibility(mi)`, median)) %>%
  ggplot(aes(weather_simple, `Visibility(mi)`, fill = weather_simple)) +
  geom_violin(trim = TRUE, scale = "width", color = NA, alpha = .8) +
  geom_boxplot(width = .2, alpha = .1, color = "grey80", outlier.shape = NA) +
  geom_point(aes(color = weather_simple),
             position = position_jitter(width = .2, seed = 1), alpha = .2, size = 1) +
  stat_summary(fun = median, geom = "point", size = 5, shape = 23,
               fill = "white", color = "black") +
  scale_fill_manual(values = pal) + scale_color_manual(values = pal) +
  coord_flip(ylim = c(0, 15)) +
  labs(title = "PHÂN BỐ TẦM NHÌN THEO TRẠNG THÁI THỜI TIẾT",
       subtitle = "Raincloud plot sắp xếp theo trung vị Tầm nhìn",
       x = NULL, y = "Tầm nhìn (miles)") +
  theme_minimal(base_size = 10) +
  theme(plot.title = element_text(face = "bold", size = 14, colour = "#4A7C59"),
        plot.subtitle = element_text(colour = "grey40"),
        legend.position = "none",
        panel.grid.major.y = element_blank(),
        panel.grid.minor = element_blank())

Câu lệnh: (1) Tạo bảng màu pastel cho từng trạng thái thời tiết. (2–4) Lọc dữ liệu. (5) Sắp xếp lại thứ tự các trạng thái thời tiết theo trung vị tầm nhìn để vẽ theo thứ tự hợp lý (từ thấp → cao). (6–7) Khởi tạo biểu đồ với tung là Tầm nhìn (mi), hoành là trạng thái thời tiết, tô màu theo nhóm thời tiết. (8) Vẽ violin plot biểu diễn mật độ phân bố tầm nhìn, không kéo dài ra ngoài phạm vi dữ liệu thật. (9) Thêm boxplot lồng bên trong để hiển thị vùng tứ phân vị. (10–11) Rải điểm dữ liệu thực để thấy mật độ; tô mờ để tránh chồng. (12–13) Đánh dấu trung vị bằng hình thoi trắng viền đen. (14–15) Áp dụng bảng màu pastel cho fill và outline. (16) Lật trục để các trạng thái thời tiết nằm dọc; giới hạn tầm nhìn từ 0–15 miles. (17–20) Thêm tiêu đề, phụ đề (Raincloud plot), và nhãn trục “Tầm nhìn (miles)”. (21–25) Dùng theme minimal, tiêu đề đậm màu xanh lá, ẩn chú giải và các đường lưới để đồ thị gọn, sáng và dễ nhìn.

Nhận xét:

Biểu đồ thể hiện sự thay đổi tầm nhìn (Visibility) dưới các điều kiện thời tiết khác nhau. Nhìn chung, tầm nhìn giảm dần từ trời quang → nhiều mây → mưa → sương mù → tuyết, phản ánh rõ tác động tiêu cực của thời tiết xấu đến khả năng quan sát khi lái xe.

Trong điều kiện trời quang, tầm nhìn đạt mức cao nhất (10–15 miles) và ổn định, thể hiện điều kiện quan sát tốt nhất. Khi nhiều mây hoặc có mưa, tầm nhìn giảm xuống còn khoảng 5–10 miles, độ phân tán tăng do ảnh hưởng của độ ẩm và mưa. Sương mù cho thấy tầm nhìn giảm mạnh chỉ còn 1–3 miles, và tuyết cũng ở mức thấp tương tự (~2–4 miles) nhưng dao động rộng hơn do phụ thuộc cường độ tuyết và quá trình tan băng.

Tóm lại, tầm nhìn bị hạn chế nghiêm trọng trong điều kiện sương mù và tuyết, trung bình khi mưa, và ổn định nhất trong thời tiết quang đãng. Biểu đồ thể hiện rõ ảnh hưởng tiêu cực của thời tiết xấu đến an toàn giao thông.

Cụm khí tượng: Độ ẩm và Áp suất

set.seed(1)
pdat <- Acc_raw %>%
  filter(between(`Pressure(in)`, 28, 31), between(`Humidity(%)`, 0, 100)) %>%
  dplyr::slice_sample(n = 200000)
ggplot(pdat, aes(`Pressure(in)`, `Humidity(%)`)) +
  geom_density_2d_filled(contour_var = "ndensity") +   # dải màu theo mật độ
  labs(title = "MẬT ĐỘ CỦA ĐỘ ẨM VÀ ÁP SUẤT",
       x = "Áp suất (inHg)", y = "Độ ẩm (%)", fill = "Mật độ") +
  theme_minimal()

Câu lệnh: (1) Cố định seed để kết quả lấy mẫu và tô mật độ tái lập. (2–4) Tạo dữ liệu pdat: lọc áp suất trong 28–31 inHg và độ ẩm 0–100%, rồi lấy mẫu 200.000 dòng để vẽ nhanh. (5) Khởi tạo biểu đồ: trục X là áp suất (inHg), trục Y là độ ẩm (%). (6) Vẽ mật độ 2D dạng vùng tô; màu thể hiện mật độ chuẩn hoá (ndensity), càng đậm → mật độ quan sát càng cao. (7–8) Thêm tiêu đề và nhãn trục, đặt thang chú giải là “Mật độ”. (9) Áp theme tối giản để biểu đồ sạch và nhẹ mắt.

Nhận xét:

Tâm mật độ cao (vùng sáng vàng–xanh nhạt) tập trung quanh áp suất ~29.8–30.0 inHg và độ ẩm 60–80%, cho thấy phần lớn quan sát nằm trong điều kiện áp suất trung bình và độ ẩm cao vừa phải.

Vùng mật độ giảm dần ra ngoài (tím đậm) thể hiện ít quan sát hơn ở áp suất thấp hoặc độ ẩm cực trị (rất khô hay ẩm bão hoà).

Cấu trúc đối xứng dọc trục áp suất gợi ý mối quan hệ yếu hoặc không tuyến tính giữa áp suất và độ ẩm — chủ yếu các thay đổi độ ẩm xảy ra trong phạm vi áp suất hẹp quanh mức trung bình khí quyển.

=> Hầu hết các hiện tượng thời tiết xảy ra quanh áp suất trung bình và độ ẩm cao vừa phải, trong khi các cực trị áp suất, độ ẩm thấp hoặc cao thường gắn liền với thời tiết bất thường hoặc nguy hiểm.

Nhiệt độ theo trạng thái thời tiết

w10  <- Acc_raw |> count(weather_simple, sort=TRUE) |> slice_head(n=10) |> pull(weather_simple)
dat  <- Acc_raw |>
  filter(weather_simple %in% w10, is.finite(temp_c)) |>
  filter(between(temp_c, quantile(temp_c,.01,na.rm=TRUE), quantile(temp_c,.99,na.rm=TRUE)))
ord  <- dat |> group_by(weather_simple) |> summarise(med = median(temp_c), .groups="drop") |> arrange(med) |> pull(weather_simple)
dat$weather_simple <- factor(dat$weather_simple, levels = ord)
meds <- dat |> group_by(weather_simple) |> summarise(med = median(temp_c), .groups="drop")
ggplot(dat, aes(temp_c, weather_simple, fill = after_stat(x))) +
  geom_density_ridges_gradient(scale = 2, rel_min_height = .003,
                               size = .25, color = "grey40",
                               quantile_lines = TRUE, quantiles = c(.05,.5,.95)) +
  geom_point(data = meds, aes(med, weather_simple), inherit.aes = FALSE,
             shape = 21, size = 2.6, fill = "white", color = "black", stroke = .5) +
  scale_fill_viridis_c(option = "inferno", name = "°C") +
  labs(title = "PHÂN BỐ NHIỆT ĐỘ THEO TRẠNG THÁI THỜI TIẾT",
       x = "Nhiệt độ (°C)", y = "Trạng thái thời tiết") +
  theme_classic() +
  theme(plot.title = element_text(face="bold"))

Câu lệnh: (1) Lấy top 10 trạng thái thời tiết xuất hiện nhiều nhất (w10). (2–4) Lọc dữ liệu: chỉ giữ bản ghi có trạng thái thuộc w10, nhiệt độ hữu hạn, và cắt đuôi 1–99% để loại ngoại lệ. (5–6) Tính trung vị nhiệt độ theo trạng thái, dùng để sắp xếp trục y (nhóm có trung vị thấp → cao). (7) Tạo bảng trung vị mỗi nhóm để vẽ điểm đánh dấu. (8) Khởi tạo đồ thị: trục x = nhiệt độ (°C), trục y = trạng thái, màu tô theo giá trị x. (9–11) Vẽ ridgeline gradient biểu diễn phân bố; kèm vạch bách phân vị 5/50/95% (đường mảnh ở giữa là median). (12–13) Thêm điểm tròn tại trung vị của từng nhóm (chấm trắng viền đen) để dễ đọc. (14) Dải màu viridis “inferno” cho thang nhiệt độ (nhãn “°C”). (15) Tiêu đề + nhãn trục. (16–18) Theme classic, làm đậm tiêu đề.

Nhận xét

Tuyết: Phân bố nhiệt tập trung ở mức rất thấp (≈ –10 → 5 °C), trung vị thấp nhất; phản ánh điều kiện lạnh sâu, điển hình mùa đông ôn đới.

Sương mù: Nhiệt độ thấp (~0–15 °C), dốc phân bố hẹp ⇒ xảy ra khi ẩm cao và không khí lạnh ổn định.

Mưa: Trải rộng từ 0–30 °C, trung vị khoảng 15–20 °C ⇒ mưa xuất hiện ở vùng ôn đến ấm, phù hợp quy luật hơi nước ngưng tụ mạnh khi nhiệt tăng.

Trời quang / Nhiều mây: Phân bố thiên phải, nhiệt độ trung bình cao hơn (20–35 °C) ⇒ thường gặp khi trời ấm, bức xạ mặt trời mạnh.

Khác: Phân bố tương tự “nhiều mây” nhưng mở rộng hơn về phía nhiệt cao ⇒ bao gồm các hiện tượng nhẹ như gió hoặc mây mỏng.

=> Biểu đồ phản ánh mối quan hệ thuận giữa nhiệt độ và trạng thái thời tiết: thời tiết càng “lạnh và ẩm” thì nhiệt độ trung bình càng thấp (tuyết, sương mù), còn thời tiết “khô – ít mây” thì nhiệt cao hơn (trời quang, khác).

Tầm nhìn và mức độ nghiêm trọng

Tỷ lệ tai nạn nặng (sev≥3) theo thời tiết và ca

tab <- within(Acc_raw, {sev3 <- Severity>=3})
heat <- aggregate(sev3 ~ weather_simple + time_period, tab, mean, na.rm=FALSE)
ggplot(heat, aes(time_period, weather_simple, fill=sev3)) +
  geom_tile() +
  scale_fill_viridis_c(labels=scales::percent) +
  labs(title="TỶ LỆ TAI NẠN NẶNG (SE≥3) THEO THỜI TIẾT VÀ CA",
       fill="Tỷ lệ") +
  theme_minimal()

Câu lệnh: (1) Tạo biến nhị phân sev3 trong bảng gốc: đánh dấu tai nạn nặng (Severity ≥ 3) = 1, còn lại = 0. (2) Gộp theo thời tiết × ca trong ngày và tính trung bình của sev3 ⇒ chính là tỷ lệ tai nạn nặng ở từng ô (vì mean của 0/1 = tỷ lệ). (3–4) Khởi tạo biểu đồ và vẽ heatmap: trục X = ca, trục Y = thời tiết, màu ô = tỷ lệ tai nạn nặng. (5) Dùng Viridis cho thang màu và hiển thị nhãn theo %. (6–7) Đặt tiêu đề và nhãn chú giải (“Tỷ lệ”). (8) Áp dụng theme tối giản cho bố cục sạch, dễ đọc.

Nhận xét:

Thời tiết: Tỷ lệ tai nạn nặng cao nhất ở mưa và “khác” (≈10–11%), trung bình ở trời quang, nhiều mây (~7–8%), thấp nhất khi sương mù (~5–6%) do xe đi chậm hơn.

Ca trong ngày: Ban đêm có tỷ lệ tai nạn nặng cao nhất, dễ hiểu vì tầm nhìn kém, mệt mỏi và thiếu ánh sáng. Buổi chiều thấp nhất, do mật độ xe cao nhưng tốc độ giảm.

Tổng thể: Tai nạn nặng chiếm khoảng 5–11%, tập trung ở điều kiện xấu hoặc ban đêm.

Mức độ nghiêm trọng theo tầm nhìn và ca trong ngày

Acc_raw %>%
  transmute(Severity, vis = pmin(`Visibility(mi)`, 50), time_period) %>%
  filter(is.finite(Severity), is.finite(vis), !is.na(time_period)) %>%
  ggplot(aes(vis, Severity)) +
  stat_binhex(bins = 35) +
  geom_smooth(method = "gam", formula = y ~ s(x, k = 6),
              se = FALSE, color = "white", linewidth = .6) +
  scale_fill_viridis_c(name = "Số vụ", trans = "sqrt") +
  facet_wrap(~ time_period, ncol = 2) +
  labs(title = "MỨC ĐỘ NGHIÊM TRỌNG THEO TẦM NHÌN TÁCH THEO CA",
       x = "Tầm nhìn (dặm)", y = "Severity") +
  theme_minimal(base_family = "Times New Roman") +
  theme(plot.title = element_text(face = "bold"))

Câu lệnh: (1–3): Lấy dữ liệu từ Acc_raw, chọn các biến Severity, Visibility(mi) (giới hạn tối đa 50 dặm để loại ngoại lệ) và time_period, đồng thời loại bỏ giá trị NA hoặc vô hạn. (4–5): Tạo biểu đồ tần suất hai chiều bằng stat_binhex(bins=35), chia không gian thành 35 ô lục giác thể hiện số vụ tai nạn theo từng mức tầm nhìn và độ nghiêm trọng. (6–7): Thêm đường xu hướng trơn (GAM) mô tả mối quan hệ trung bình giữa tầm nhìn và mức độ nghiêm trọng, bỏ dải sai số (se=FALSE), màu trắng để nổi bật. (8–9): Dùng thang màu Viridis để biểu diễn mật độ (chuyển đổi căn bậc hai) và facet_wrap chia biểu đồ theo từng ca trong ngày (2 cột). (10–13): Gán tiêu đề, nhãn trục và định dạng giao diện Times New Roman, tiêu đề in đậm cho đẹp và rõ ràng.

Nhận xét

Phân bố chung: Phần lớn các vụ tai nạn tập trung ở vùng tầm nhìn ngắn (0–10 dặm), cho thấy tầm nhìn kém là yếu tố rủi ro lớn.

Theo mức độ nghiêm trọng (Severity): Các vụ mức 2 (trung bình) chiếm tỷ lệ cao nhất, trong khi mức 3–4 (nặng) xuất hiện thưa hơn, chủ yếu khi tầm nhìn rất thấp.

Theo ca: Chiều và Sáng có mật độ tai nạn cao nhất, có thể do lưu lượng giao thông lớn (giờ đi làm và tan ca). Đêm và Tối ít vụ hơn nhưng có xu hướng mức độ nghiêm trọng cao hơn — gợi ý ảnh hưởng của thiếu sáng và mệt mỏi khi lái xe.

Không gian và cấu trúc liên hệ

Mức độ nghiêm trọng với quãng đường theo bang (Top 8)

s8 <- names(sort(table(Acc_raw$State),decreasing=TRUE))[1:8]
ggplot(subset(Acc_raw, State %in% s8),
       aes(`Distance(mi)`, Severity)) +
  geom_point(alpha=.06) +
  geom_smooth(se=FALSE) +
  facet_wrap(~State) +
  geom_hline(yintercept=median(Acc_raw$Severity, na.rm=TRUE),
             linetype="dashed") +
  labs(title="MỨC ĐỘ NGHIÊM TRỌNG VÀ QUÃNG ĐƯỜNG SỰ CỐ") +
  theme_minimal()

Câu lệnh: (1) Lấy 8 bang có nhiều tai nạn nhất bằng cách đếm tần suất, sắp xếp giảm dần và chọn 8 giá trị đầu. (2–3) Tạo biểu đồ từ dữ liệu chỉ gồm 8 bang top đầu, trục X là Distance (mi) (độ dài đoạn đường bị ảnh hưởng), trục Y là Severity (mức độ nghiêm trọng). (4–5) Vẽ các điểm dữ liệu phân tán và đường xu hướng trơn thể hiện mối quan hệ giữa quãng đường và mức độ nghiêm trọng. (6) Tách biểu đồ riêng cho từng bang. (7) Thêm đường ngang tại trung vị toàn quốc của Severity để làm mốc tham chiếu (nét gạch dashed). (8–10) Đặt tiêu đề và dùng giao diện tối giản.

Nhận xét:

Biểu đồ mô tả mối quan hệ giữa độ nghiêm trọng (Severity) và quãng đường sự cố (Distance) ở tám bang có số liệu nhiều nhất. Các điểm đen biểu thị dữ liệu thực tế, đường xanh thể hiện xu hướng trung bình (loess) và đường đen nét đứt cho biết mức nghiêm trọng trung bình khoảng 2.

Phần lớn tai nạn tập trung ở mức Severity = 2 và Distance < 20 dặm, phản ánh rằng đa số là sự cố nhẹ trong đô thị. Ở hầu hết các bang, đường xu hướng dao động quanh mức 2, cho thấy quãng đường sự cố không ảnh hưởng đáng kể đến mức độ nghiêm trọng.

Một số bang như Texas (TX) và Florida (FL) có xu hướng đường xanh tăng nhẹ khi quãng đường dài hơn, hàm ý rằng các vụ tai nạn xảy ra trên đường cao tốc hoặc tuyến xa lộ thường nghiêm trọng hơn. Ngược lại, các bang còn lại duy trì xu hướng ổn định, cho thấy đặc trưng giao thông đô thị chiếm ưu thế.

Ma trận tương quan

df_num <- Acc_raw %>%
  transmute(
    Severity = as.numeric(Severity),
    Distance_mi = as.numeric(`Distance(mi)`),
    TempC = coalesce(as.numeric(temp_c), (as.numeric(`Temperature(F)`) - 32) * 5/9),
    Humidity = as.numeric(`Humidity(%)`),
    Pressure = as.numeric(`Pressure(in)`),
    Visibility = as.numeric(`Visibility(mi)`)
  ) %>%
  drop_na()
C <- cor(df_num, use = "pairwise.complete.obs")
ord <- hclust(as.dist(1 - abs(C)))$order
C2 <- C[ord, ord]
cor_long <- as.data.frame(C2) %>%
  tibble::rownames_to_column("Var1") %>%
  pivot_longer(-Var1, names_to = "Var2", values_to = "r") %>%
  mutate(
    Var1 = factor(Var1, levels = rownames(C2)),
    Var2 = factor(Var2, levels = colnames(C2))) %>%
  filter(as.integer(Var1) <= as.integer(Var2)) %>%
  mutate(r_lab = scales::label_number(accuracy = 0.01, decimal.mark = ",")(r))
ggplot(cor_long, aes(Var2, Var1, fill = r)) +
  geom_tile(color = "grey90", linewidth = 0.3) +
  geom_text(aes(label = r_lab), size = 3) +
  scale_fill_gradient2(limits = c(-1, 1), midpoint = 0,
                       low = "#2b8cbe", mid = "white", high = "#de2d26", name = "r") +
  coord_fixed() +
  labs(title = "MA TRẬN TƯƠNG QUAN") +
  theme_minimal(base_size = 12) +
  theme(
    panel.grid = element_blank(),
    axis.text.x = element_text(angle = 45, hjust = 1),
    legend.position = "right")

Câu lệnh: (1–9) Chuyển các cột định lượng sang dạng số gồm Severity, Distance(mi), Temperature, Humidity, Pressure, Visibility, trong đó nhiệt độ đổi từ °F sang °C để thống nhất đơn vị. Loại bỏ các hàng bị thiếu dữ liệu để thu được bảng df_num sạch. (10–12) Tính ma trận tương quan Pearson giữa các biến định lượng bằng cor(), sau đó gom nhóm các biến có mối liên hệ tương tự bằng phân cụm thứ bậc (hclust) và sắp xếp lại thứ tự biến (C2 <- C[ord, ord]) để biểu đồ gọn và dễ quan sát hơn. (13–20) Chuyển ma trận C2 sang dạng dài (long format), trong đó mỗi dòng biểu diễn một cặp biến và giá trị hệ số r. Giữ lại một nửa ma trận (tam giác trên) để tránh trùng lặp, đồng thời tạo nhãn r_lab làm tròn giá trị tương quan đến 0,01. (21–27) Vẽ heatmap biểu diễn tương quan: geom_tile() tạo các ô màu, geom_text() hiển thị hệ số r. Thang màu được đặt với xanh cho âm mạnh, trắng cho trung tính (0), và đỏ cho dương mạnh. (28–33) Cố định tỷ lệ ô vuông , thêm tiêu đề và nhãn trục. Áp dụng theme_minimal() để làm gọn giao diện, ẩn đường lưới, xoay nhãn trục X 45°, và đặt chú giải ở bên phải biểu đồ.

Nhận xét:

Quan hệ giữa các biến nhìn chung yếu, hầu hết các hệ số r đều nằm quanh 0 → cho thấy các yếu tố không tuyến tính mạnh với nhau.

Nhiệt độ (TempC) có mối tương quan dương nhẹ với Độ ẩm (r = 0,30) và âm với Áp suất (r = –0,33), phù hợp về mặt khí tượng: khi nhiệt độ tăng, áp suất thường giảm và độ ẩm thay đổi cùng chiều.

Độ ẩm và Tầm nhìn (r = –0,40) có tương quan âm trung bình, phản ánh hiện tượng thực tế: độ ẩm cao (mưa, sương mù) làm giảm tầm nhìn.

Severity (mức độ nghiêm trọng) gần như không tương quan đáng kể với các yếu tố khí tượng hay khoảng cách (|r| < 0,1), hàm ý rằng độ nghiêm trọng của tai nạn phụ thuộc nhiều hơn vào yếu tố hành vi, mật độ giao thông hoặc điều kiện đường hơn là thời tiết riêng lẻ.

Áp suất và các biến khác chỉ có tương quan yếu, chủ yếu quanh ±0,1 → tác động không đáng kể trong dữ liệu.

Bộ dữ liệu được sử dụng trong nghiên cứu này được thu thập từ báo cáo tài chính hợp nhất của Ngân hàng TMCP Sài Gòn – Hà Nội (SHB) trong giai đoạn 2014–2023. Dữ liệu bao gồm các chỉ tiêu trọng yếu phản ánh quy mô tài sản, cơ cấu nguồn vốn, hiệu quả hoạt động và khả năng sinh lời, phục vụ cho việc đánh giá mức độ tăng trưởng và ổn định tài chính của ngân hàng qua thời gian. Bộ dữ liệu gồm 13 biến đại diện cho các nhóm thông tin chính:

Nhóm đặc điểm thời gian: biến Năm thể hiện giai đoạn quan sát, giúp theo dõi sự biến động qua các kỳ.

Nhóm quy mô và cấu trúc tài sản – nguồn vốn: bao gồm Tổng cộng tài sản, Vốn chủ sở hữu, Nợ phải trả, Tiền gửi tại NHNN, Tiền & vàng gửi tại các TCTD khác, Chứng khoán kinh doanh, và Tiền mặt & chứng từ có giá trị. Các biến này phản ánh năng lực tài chính, mức độ an toàn vốn và chính sách sử dụng tài sản của SHB.

Nhóm kết quả kinh doanh: gồm Doanh thu thuần, Thu nhập lãi, Thu nhập lãi thuần, Thu nhập ngoài lãi, và Lợi nhuận sau thuế. Đây là các biến phản ánh trực tiếp hiệu quả hoạt động, khả năng tạo lợi nhuận và năng lực sinh lời của ngân hàng.

Tất cả dữ liệu được xử lý và tổng hợp theo đơn vị tiền tệ là đồng Việt Nam (VND), dưới dạng giá trị tuyệt đối, được trích xuất từ các báo cáo tài chính đã kiểm toán của ngân hàng. Quá trình thu thập và xử lý dữ liệu được thực hiện thủ công kết hợp với công cụ R nhằm đảm bảo tính chính xác và đồng nhất giữa các năm.

Bộ dữ liệu tài chính của SHB đóng vai trò là nguồn thông tin định lượng cốt lõi giúp nhóm nghiên cứu phân tích xu hướng tăng trưởng và mức độ ổn định tài chính, từ đó rút ra những nhận định khoa học và có giá trị thực tiễn đối với hoạt động ngân hàng thương mại Việt Nam trong giai đoạn nghiên cứu.

Giới thiệu bộ dữ liệu

Giới thiệu và tổng quan bộ dữ liệu

Đọc bộ dữ liệu

ds <- read_csv("~/Library/CloudStorage/OneDrive-UFM/P2bctc.csv")

Câu lệnh: (1) câu lệnh thuộc gói readr được dùng để đọc tệp dữ liệu định dạng CSV.

Xem kích thước bộ dữ liệu

dim(ds)

## [1] 10 13

Câu lệnh: (1) trả về một vector gồm 2 số: số dòng (quan sát) và số cột (biến) của đối tượng dữ liệu ds.

Kết quả: [1] 10 13 nghĩa là bảng dữ liệu có 10 dòng và 13 cột → 10 quan sát, 13 biến.

Hiển thị 5 dòng đầu tiên của dữ liệu

ds_h <- head(ds, 5) |>
  mutate(across(where(is.numeric) & !matches("^Năm$"),
                ~ prettyNum(.x, big.mark=".", decimal.mark=",",
                            preserve.width="none"))) |>
  as.data.frame()
ds_h

Câu lệnh: (1) Tạo bản xem nhanh gồm 5 dòng đầu của dữ liệu. (2–4) Định dạng mọi cột số (trừ cột tên “Năm”): thêm dấu chấm ngăn cách hàng nghìn và dấu phẩy cho thập phân; không giữ độ rộng ký tự. (5) Đổi về data.frame để in hiển thị đơn giản. (6) In đối tượng đã xử lý.

Nhận xét:

Bảng trên cho thấy dữ liệu đã được định dạng lại dễ đọc, với dấu chấm ngăn cách hàng nghìn. Các giá trị tăng dần theo năm, thể hiện quy mô tài sản, vốn và lợi nhuận của SHB liên tục mở rộng. Điều này phản ánh xu hướng tăng trưởng ổn định của ngân hàng trong giai đoạn 2019–2023.

Lược đồ và kiểu dữ liệu

Kiểm tra cấu trúc và kiểu dữ liệu của bộ dữ liệu

glimpse(ds)

## Rows: 10
## Columns: 13
## $ Năm                                 <dbl> 2023, 2022, 2021, 2020, …
## $ `Tổng cộng tài sản`                 <dbl> 630500685000000, 5509041…
## $ `Vốn chủ sở hữu`                    <dbl> 50098280000000, 42904471…
## $ `Lợi nhuận sau thuế`                <dbl> 7324758000000, 772891800…
## $ `Thu nhập lãi`                      <dbl> 58898319000000, 31743193…
## $ `Tiền mặt & Chứng từ có giá trị`    <dbl> 1370849000000, 189754500…
## $ `Tiền gửi tại NHNN`                 <dbl> 54763646000000, 15145862…
## $ `Tiền & vàng gửi tại các TCTD khác` <dbl> 63548928000000, 63008862…
## $ `Chứng khoán kinh doanh`            <dbl> 7792742000000, 154700000…
## $ `Thu Nhập Ngoài Lãi`                <dbl> 1299873000000, 178667200…
## $ `Nợ phải trả`                       <dbl> 580402405000000, 4995236…
## $ `Thu nhập lãi thuần`                <dbl> 19285353000000, 17550084…
## $ `Doanh thu thuần`                   <dbl> 21328397000000, 19340982…

Câu lệnh: (1) In tóm tắt cấu trúc: số dòng/cột, rồi liệt kê tên biến → kiểu dữ liệu (ví dụ ) kèm mẫu vài giá trị đầu cho mỗi biến. Hữu ích để nhìn nhanh sơ lược bộ dữ liệu.

Kết quả

Bộ dữ liệu có 10 dòng và 13 cột.

Tất cả biến đang là số thực : Năm, Tổng cộng tài sản, Vốn chủ sở hữu, Lợi nhuận sau thuế, Thu nhập lãi, Tiền mặt, Tiền gửi tại NHNN, Tiền, vàng gửi, Chứng khoán kinh doanh, Thu Nhập Ngoài Lãi, Nợ phải trả, Thu nhập lãi thuần, Doanh thu thuần.

Cột Năm hiển thị dãy giảm dần (2023 →2014) và đơn vị là đồng.

Kiểm tra Số lượng cột dạng số với không phải số

nums  <- ds %>% dplyr::select_if(is.numeric) %>% ncol()
c(numeric = nums, non_numeric = ncol(ds) - nums)

##     numeric non_numeric 
##          13           0

Câu lệnh: (1) Đếm số cột kiểu số trong bảng ds. (2) Tạo vector gồm hai giá trị: số cột numeric và số cột không phải numeric.

Kết quả: Trả về thống kê nhanh 13 cột numeric và non-numeric của ds.

Kiểm tra kiểu dữ liệu

sapply(ds, class)

##                               Năm                 Tổng cộng tài sản 
##                         "numeric"                         "numeric" 
##                    Vốn chủ sở hữu                Lợi nhuận sau thuế 
##                         "numeric"                         "numeric" 
##                      Thu nhập lãi    Tiền mặt & Chứng từ có giá trị 
##                         "numeric"                         "numeric" 
##                 Tiền gửi tại NHNN Tiền & vàng gửi tại các TCTD khác 
##                         "numeric"                         "numeric" 
##            Chứng khoán kinh doanh                Thu Nhập Ngoài Lãi 
##                         "numeric"                         "numeric" 
##                       Nợ phải trả                Thu nhập lãi thuần 
##                         "numeric"                         "numeric" 
##                   Doanh thu thuần 
##                         "numeric"

Câu lệnh: (1) Áp dụng hàm class lên mỗi cột của ds và trả về vector các lớp dữ liệu tương ứng (ở đây in theo tên cột).

Kết quả: Kết quả cho thấy tất cả các biến trong bộ dữ liệu đều thuộc kiểu số (numeric). Điều này phù hợp vì các cột đều là chỉ tiêu tài chính như tài sản, vốn, lợi nhuận, thu nhập… cần định dạng số để tính toán, so sánh và trực quan hóa trong các phân tích tiếp theo.

Tên biến và tính hợp lệ

Liệt kê tên các biến

names(ds)

##  [1] "Năm"                              
##  [2] "Tổng cộng tài sản"                
##  [3] "Vốn chủ sở hữu"                   
##  [4] "Lợi nhuận sau thuế"               
##  [5] "Thu nhập lãi"                     
##  [6] "Tiền mặt & Chứng từ có giá trị"   
##  [7] "Tiền gửi tại NHNN"                
##  [8] "Tiền & vàng gửi tại các TCTD khác"
##  [9] "Chứng khoán kinh doanh"           
## [10] "Thu Nhập Ngoài Lãi"               
## [11] "Nợ phải trả"                      
## [12] "Thu nhập lãi thuần"               
## [13] "Doanh thu thuần"

Câu lệnh: (1) Trả về vector tên các cột (biến) của bảng dữ liệu ds.

Kết quả:

Năm (Year): Biến định danh thời gian, thể hiện năm tài chính tương ứng với các giá trị số liệu. Đây là biến cơ sở giúp sắp xếp và theo dõi xu hướng biến động của các chỉ tiêu tài chính qua từng giai đoạn.

Tổng cộng tài sản (Total Assets): Phản ánh toàn bộ quy mô tài sản mà ngân hàng sở hữu hoặc quản lý tại thời điểm cuối kỳ, bao gồm cả tài sản ngắn hạn và dài hạn. Đây là chỉ tiêu tổng hợp quan trọng để đánh giá quy mô, năng lực tài chính và mức độ mở rộng hoạt động của ngân hàng.

Vốn chủ sở hữu (Owner’s Equity): Đại diện cho phần vốn thuộc quyền sở hữu của cổ đông và các quỹ nội bộ sau khi trừ đi nợ phải trả. Chỉ tiêu này thể hiện năng lực tự chủ tài chính và khả năng chống chịu rủi ro của ngân hàng.

Lợi nhuận sau thuế (Net Profit After Tax):Là phần lợi nhuận ròng còn lại sau khi đã trừ toàn bộ chi phí và thuế thu nhập doanh nghiệp. Đây là biến phản ánh kết quả cuối cùng của hoạt động kinh doanh, đồng thời được dùng để tính các chỉ tiêu sinh lời như ROA, ROE.

Thu nhập lãi (Interest Income):Biến này thể hiện tổng thu nhập từ các hoạt động cho vay và đầu tư sinh lãi – nguồn thu chủ yếu của ngân hàng thương mại. Mức độ tăng trưởng của chỉ tiêu này phản ánh hiệu quả khai thác nguồn vốn và khả năng mở rộng tín dụng.

Tiền mặt và chứng từ có giá trị (Cash and Valuable Papers):Phản ánh lượng tiền mặt, vàng bạc, ngoại tệ, cùng các loại chứng từ có giá trị như séc, kỳ phiếu… mà ngân hàng đang nắm giữ. Đây là tài sản có tính thanh khoản cao, giúp đảm bảo khả năng chi trả trong ngắn hạn.

Tiền gửi tại NHNN (Deposits at State Bank):Là khoản tiền ngân hàng gửi tại Ngân hàng Nhà nước nhằm đáp ứng yêu cầu dự trữ bắt buộc, thanh toán liên ngân hàng hoặc các nghĩa vụ khác. Chỉ tiêu này phản ánh mối quan hệ và nghĩa vụ thanh khoản của ngân hàng với cơ quan quản lý tiền tệ.

Tiền và vàng gửi tại các TCTD khác (Deposits at Other Credit Institutions):Thể hiện khoản tiền, vàng hoặc ngoại tệ ngân hàng gửi tại các tổ chức tín dụng khác nhằm mục đích thanh toán, dự phòng hoặc đầu tư ngắn hạn. Chỉ tiêu này cho thấy mức độ liên kết và hợp tác trong hệ thống ngân hàng.

Chứng khoán kinh doanh (Trading Securities):Là giá trị các khoản đầu tư chứng khoán mà ngân hàng nắm giữ nhằm mục tiêu mua bán ngắn hạn kiếm lời. Biến này phản ánh khả năng đầu tư và linh hoạt trong chiến lược sử dụng vốn.

Thu nhập ngoài lãi (Non-interest Income):Biểu thị nguồn thu từ các hoạt động phi tín dụng như dịch vụ thanh toán, kinh doanh ngoại hối, bảo hiểm, tư vấn tài chính… Chỉ tiêu này cho thấy mức độ đa dạng hóa nguồn thu của ngân hàng, giúp giảm phụ thuộc vào hoạt động tín dụng.

Nợ phải trả (Liabilities):Phản ánh tổng giá trị nghĩa vụ tài chính mà ngân hàng phải thanh toán cho các bên liên quan, bao gồm tiền gửi của khách hàng, vay liên ngân hàng và các khoản nợ khác. Đây là chỉ tiêu quan trọng trong đánh giá cơ cấu vốn và rủi ro tài chính.

Thu nhập lãi thuần (Net Interest Income):Là phần chênh lệch giữa thu nhập lãi và chi phí lãi, thể hiện hiệu quả hoạt động cốt lõi trong kinh doanh tín dụng. Mức tăng của biến này cho thấy khả năng sinh lợi từ hoạt động cho vay và đầu tư tài chính.

Doanh thu thuần (Net Revenue):Phản ánh tổng doanh thu sau khi đã trừ các khoản giảm trừ, chiết khấu và hoàn nhập. Đây là chỉ tiêu tổng hợp quan trọng để đo lường quy mô hoạt động kinh doanh và hiệu quả bán hàng của ngân hàng qua từng năm.

Kiểm tra tên cột trùng/trống

list(ten_trung = names(ds)[duplicated(names(ds))],
     co_ten_trong = any(names(ds)==""), 
     so_cot = ncol(ds))

## $ten_trung
## character(0)
## 
## $co_ten_trong
## [1] FALSE
## 
## $so_cot
## [1] 13

Câu lệnh: (1) Liệt kê các tên cột bị trùng trong bảng. (2) Trả về TRUE/FALSE: có cột tên rỗng (““) hay không. (3) Đếm số cột của bảng.

Kết quả: Không có tên cột nào trùng và trống. Tổng số cột trong bộ dữ liệu là 13 cột gồm 13 biến.

Chất lượng dữ liệu

Kiểm tra NA cho mỗi cột

colSums(is.na(ds))

##                               Năm                 Tổng cộng tài sản 
##                                 0                                 0 
##                    Vốn chủ sở hữu                Lợi nhuận sau thuế 
##                                 0                                 0 
##                      Thu nhập lãi    Tiền mặt & Chứng từ có giá trị 
##                                 0                                 0 
##                 Tiền gửi tại NHNN Tiền & vàng gửi tại các TCTD khác 
##                                 0                                 0 
##            Chứng khoán kinh doanh                Thu Nhập Ngoài Lãi 
##                                 0                                 0 
##                       Nợ phải trả                Thu nhập lãi thuần 
##                                 0                                 0 
##                   Doanh thu thuần 
##                                 0

Câu lệnh: (1) Lệnh trên đếm tổng số giá trị thiếu (NA) ở từng cột của bảng Accidents

Kết quả: Tất cả các cột đều trả về 0 ⇒ không có NA trong toàn bộ 13 biến. Dữ liệu hoàn chỉnh, không cần bước xử lý thiếu (impute/loại dòng).

Kiểm tra giá trị trùng lặp

sum(duplicated(ds))

## [1] 0

Câu lệnh: (1) sẽ kiểm tra các dòng trùng lặp trong toàn bộ bảng dữ liệu ds.

Kết quả: Kết quả 0, không có dòng trùng lặp.

Xử lý dữ liệu thô và mã hoá dữ liệu

Sao lưu bản gốc trước khi xử lý: mọi bước làm sạch/biến đổi sẽ thực hiện trên ds1, còn ds được giữ nguyên để đối chiếu hoặc khôi phục khi cần.

ds1 <- ds
ds1 <- ds1 %>% arrange(`Năm`)

Chuẩn hoá cấu trúc biến

Chuẩn kiểu biến “Năm” thành số nguyên

Để chuẩn hóa kiểu dữ liệu của biến ““Năm”” từ dạng số thực thành số nguyên, giúp dễ dàng sắp xếp, nhóm, tính toán và vẽ biểu đồ theo năm. Tránh tình trạng Năm đang là character hoặc double (ví dụ 2020.0), gây lỗi khi join/merge hoặc khi dùng.

ds1 <- ds1 %>% mutate(`Năm` = as.integer(`Năm`))

Câu lệnh: (1) Dòng lệnh này chuyển cột “Năm” trong bảng ds1 sang kiểu số nguyên (integer).

Chuẩn hoá tên biến

Các tên biến có dấu, khoảng trắng hoặc ký tự đặc biệt thường gây lỗi khi gọi trong các hàm của dplyr hoặc ggplot2. Do đó, việc chuyển toàn bộ tên biến về dạng không dấu, chữ thường, và ngăn cách bằng dấu gạch dưới (_) giúp cho việc lập trình, trực quan hoá và xuất báo cáo được mạch lạc và dễ thao tác hơn.

simple_names <- function(x){
  x <- stringi::stri_trans_general(x, "Latin-ASCII")  
  x <- tolower(x)
  x <- gsub("[^a-z0-9]+","_", x)
  gsub("(^_|_$)", "", x)}
names(ds1) <- simple_names(names(ds1))
names(ds1)

##  [1] "nam"                            
##  [2] "tong_cong_tai_san"              
##  [3] "von_chu_so_huu"                 
##  [4] "loi_nhuan_sau_thue"             
##  [5] "thu_nhap_lai"                   
##  [6] "tien_mat_chung_tu_co_gia_tri"   
##  [7] "tien_gui_tai_nhnn"              
##  [8] "tien_vang_gui_tai_cac_tctd_khac"
##  [9] "chung_khoan_kinh_doanh"         
## [10] "thu_nhap_ngoai_lai"             
## [11] "no_phai_tra"                    
## [12] "thu_nhap_lai_thuan"             
## [13] "doanh_thu_thuan"

Câu lệnh: (1) Tạo hàm chuẩn hoá tên biến. (2–5) Chuẩn hoá tên: bỏ dấu → ASCII, chuyển thường, thay ký tự lạ bằng “”, xoá “” ở đầu/cuối. (6) Áp dụng cho toàn bộ tên cột của ds1. (7) Xem lại danh sách tên mới.

Ép numeric an toàn cho các cột số (bỏ dấu phẩy nếu có)

Khi nhập từ Excel/CSV, nhiều cột số bị lưu dưới dạng ký tự vì có dấu phân tách hàng nghìn (,). Nếu không ép lại về numeric, các phép tính (tăng trưởng, thống kê mô tả, vẽ biểu đồ, hồi quy) sẽ lỗi hoặc cho kết quả sai.

num_cols <- setdiff(names(ds1), "nam")
num_cols <- intersect(num_cols, names(ds1))
ds1[num_cols] <- lapply(ds1[num_cols], \(x) as.numeric(gsub(",","",x)))

Câu lệnh: (1) Lấy tất cả tên cột trừ “nam”. (2) So khớp lại với tên cột hiện có. (3) Với các cột ấy: xóa dấu phẩy trong chuỗi số rồi ép sang numeric.

Xây dựng chỉ tiêu tài chính

Tạo biến log_total_assets

Trong nghiên cứu tài chính, biến “Tổng cộng tài sản” thường có giá trị rất lớn và chênh lệch mạnh giữa các năm, khiến dữ liệu bị lệch phải và dễ gây sai lệch khi phân tích. Việc logarit hóa tổng tài sản giúp thu hẹp độ phân tán, giảm ảnh hưởng ngoại lai và ổn định phương sai. Đồng thời, khi các biến được log hoá, hệ số hồi quy có thể diễn giải theo tỷ lệ phần trăm, giúp kết quả mang ý nghĩa kinh tế rõ ràng hơn.

ds1 <- ds1 %>% mutate(log_total_assets = log(tong_cong_tai_san))

Câu lệnh: (1) Dòng lệnh này tạo thêm biến mới tên là log_total_assets, chứa giá trị logarit tự nhiên của cột “tong_cong_tai_san”.

Kết quả là mỗi giá trị trong cột “Tổng cộng tài sản” được chuyển đổi sang giá trị logarit tương ứng, phản ánh quy mô tài sản theo thang logarit.

Tạo biến ROA, ROE

ROA (Return on Assets) và ROE (Return on Equity) là hai chỉ tiêu cốt lõi để đánh giá hiệu quả sử dụng tài sản và vốn chủ sở hữu của ngân hàng. ROA đo lường khả năng sinh lợi trên mỗi đồng tài sản, phản ánh hiệu quả quản lý tài sản và cấu trúc bảng cân đối. ROE đo lường mức sinh lời cho cổ đông, chịu tác động của đòn bẩy tài chính (vốn chủ sở hữu thấp với lợi nhuận giữ nguyên sẽ làm ROE cao). Hai chỉ tiêu này sẽ được sử dụng ở các phần sau để mô tả, so sánh theo thời gian và là biến phụ thuộc/giải thích trong các phân tích sâu hơn.

ds1 <- ds1 %>%
  mutate(
    roa = loi_nhuan_sau_thue / tong_cong_tai_san,
    roe = loi_nhuan_sau_thue / von_chu_so_huu)

Câu lệnh: (1-3) Thêm cột mức sinh lời trên tài sản. (1-4) Thêm cột mức sinh lời trên vốn CSH.

Tạo biến debt_ratio

Đánh giá mức độ đòn bẩy tài chính qua hệ số nợ (debt_ratio). Tỷ lệ cao cho thấy phụ thuộc nhiều vào vốn vay, lợi nhuận có thể cao nhưng rủi ro lớn; tỷ lệ thấp phản ánh cơ cấu vốn an toàn hơn, song lợi nhuận tiềm năng thấp hơn.

ds1 <- ds1 %>% mutate(debt_ratio = no_phai_tra / tong_cong_tai_san)

Câu lệnh: (1) Tạo thêm một biến debt_ratio mới trong bảng ds1.

Tạo biến cash_ratio

Đánh giá khả năng thanh khoản tức thời của ngân hàng. Cash_ratio cao thể hiện an toàn thanh khoản nhưng giảm hiệu quả sinh lời, trong khi tỷ lệ thấp giúp tăng lợi nhuận song tiềm ẩn rủi ro thiếu hụt tiền mặt.

ds1 <- ds1 %>% mutate(cash_ratio = tien_mat_chung_tu_co_gia_tri / tong_cong_tai_san)

Câu lệnh: (1) Tạo thêm một biến cash_ratio mới trong bảng ds1.

Tạo biến growth_assets

Đo tốc độ tăng trưởng tổng tài sản để xác định mức mở rộng hay thu hẹp quy mô ngân hàng qua các năm. Growth_assets dương cho thấy tài sản tăng, âm phản ánh thu hẹp; tăng trưởng cao kéo theo rủi ro tín dụng và thanh khoản nên cần so sánh với chỉ tiêu an toàn vốn để đánh giá bền vững.

ds1 <- ds1 %>% mutate(growth_assets =
                        (tong_cong_tai_san - lag(tong_cong_tai_san)) / 
                        lag(tong_cong_tai_san))

Câu lệnh: (1) Tạo thêm một biến growth_assets mới trong bảng ds1.

Tạo biến non_interest_income_share

Tỷ trọng thu nhập ngoài lãi thể hiện mức độ đa dạng hóa nguồn thu của ngân hàng. Giá trị cao cho thấy cơ cấu doanh thu ổn định, ít phụ thuộc vào tín dụng; ngược lại, giá trị thấp phản ánh sự phụ thuộc lớn vào chênh lệch lãi suất.

ds1 <- ds1 %>% mutate(non_interest_income_share = thu_nhap_ngoai_lai / (thu_nhap_lai + thu_nhap_ngoai_lai))

Câu lệnh: (1) Tạo thêm một biến non_interest_income_share mới trong bảng ds1.

Phân nhóm đặc trưng

Phân tổ theo quy mô tài sản

Việc phân tổ dữ liệu theo quy mô tài sản để theo dõi quá trình mở rộng quy mô và so sánh hiệu quả hoạt động giữa các giai đoạn phát triển của ngân hàng.

ds1 <- ds1 %>%
  mutate(asset_group = case_when(
    log_total_assets <= quantile(log_total_assets, 0.33, na.rm = TRUE) ~ "Nhỏ",
    log_total_assets <= quantile(log_total_assets, 0.67, na.rm = TRUE) ~ "Trung bình",
    TRUE ~ "Lớn"))

Câu lệnh: (1) Phân nhóm quy mô tài sản dựa theo log của tổng tài sản. (2–3) Xác định hai ngưỡng tại phân vị 33% và 67% để chia thành ba mức: 33% là “Nhỏ”, 33–67% là “Trung bình”, >67% là “Lớn”. (4) Các giá trị khuyết được bỏ qua khi tính phân vị.

Phân tổ theo giai đoạn và đòn bẩy

Phân chia dữ liệu theo giai đoạn và mức đòn bẩy để so sánh biến động tài chính của SHB trước – sau 2018. Cách phân tổ này giúp đánh giá mối liên hệ giữa đòn bẩy, sinh lời (ROA, ROE) và thanh khoản (cash_ratio), qua đó nhận diện mức độ rủi ro và hiệu quả trong từng giai đoạn.

ds1 <- ds1 %>% mutate(period = ifelse(nam <= 2018, "2014–2018", "2019–2023"),
         leverage_group = case_when(
           debt_ratio <= quantile(debt_ratio, 0.33, na.rm=TRUE) ~ "Thấp",
           debt_ratio <= quantile(debt_ratio, 0.67, na.rm=TRUE) ~ "Trung bình",
           TRUE ~ "Cao"))

Câu lệnh: (1) Tạo biến period: năm ≤ 2018 → “2014–2018”, ngược lại → “2019–2023”. (2) Bắt đầu biến leverage_group theo ngưỡng phân vị của debt_ratio. (3) Nếu ≤ phân vị 33% → gán “Thấp” (bỏ qua NA). (4) Nếu ≤ phân vị 67% → gán “Trung bình”. (5) Còn lại → “Cao”.

Các thống kê cơ bản về bộ dữ liệu

Thống kê mô tả cơ bản

Top 5 năm có Tổng cộng tài sản lớn nhất

ds1 %>%
  dplyr::arrange(dplyr::desc(tong_cong_tai_san)) %>%
  dplyr::select(nam, tong_cong_tai_san) %>%
  head(5) %>%
  dplyr::mutate(
    tong_cong_tai_san = scales::number(
      tong_cong_tai_san, big.mark = ".", decimal.mark = ",", accuracy = 1))

Câu lệnh:(1) Lấy dữ liệu gốc và bắt đầu chuỗi xử lý. (2) Sắp xếp các năm theo tổng tài sản giảm dần để tìm ra những năm có giá trị cao nhất. (3) Giữ lại hai biến cần thiết là năm và tổng tài sản. (4) Chỉ lấy năm dòng đầu tiên để chọn ra top 5 năm lớn nhất. (5) Định dạng lại giá trị tổng tài sản bằng cách thêm dấu chấm cho hàng nghìn, dấu phẩy cho thập phân và làm tròn đến đơn vị, giúp bảng kết quả hiển thị rõ ràng, dễ đọc hơn.

Kết quả: Kết quả cho thấy tổng tài sản của SHB tăng liên tục qua các năm, thể hiện xu hướng mở rộng quy mô rõ rệt. Cụ thể, từ năm 2019 đến 2023, tổng tài sản tăng mạnh từ khoảng 365 nghìn tỷ lên hơn 630 nghìn tỷ đồng, tương đương mức tăng gần 1,7 lần sau 5 năm. Mức tăng rõ rệt nhất nằm ở giai đoạn 2020–2023, phản ánh chiến lược tăng trưởng mạnh về quy mô và mở rộng hoạt động kinh doanh của ngân hàng.

Thống kê mô tả cho biến ROA

ds1 %>% summarise(
    mean_roa = mean(roa, na.rm = TRUE),
    sd_roa   = sd(roa, na.rm = TRUE),
    min_roa  = min(roa, na.rm = TRUE),
    q1_roa   = quantile(roa, 0.25, na.rm = TRUE),
    median_roa = median(roa, na.rm = TRUE),
    q3_roa   = quantile(roa, 0.75, na.rm = TRUE),
    max_roa  = max(roa, na.rm = TRUE))

Câu lệnh: (1-8) Tạo bảng tóm tắt cho roa về trung bình, độ lệch chuẩn, giá trị nhỏ nhất, phân vị 25%, trung vị, phân vị 75% và giá trị lớn nhất.

Kết quả: ROA trung bình của SHB giai đoạn 2014–2023 đạt 0,0071 (≈ 0,71%), nghĩa là mỗi 100 đồng tài sản tạo ra khoảng 0,71 đồng lợi nhuận sau thuế. Độ lệch chuẩn 0,0035 cho thấy biến động thấp, phản ánh hiệu quả sử dụng tài sản ổn định. Giá trị ROA dao động từ 0,0039 đến 0,0140, biên độ không lớn nhưng vẫn thể hiện xu hướng cải thiện, khi mức cao nhất gấp gần 3,6 lần mức thấp nhất. Các tứ phân vị lần lượt là 0,0048 – 0,0058 – 0,0091, cho thấy phần lớn các năm ngân hàng duy trì ROA trong khoảng 0,48–0,91%, nằm trong vùng an toàn và hiệu quả ổn định của hệ thống ngân hàng Việt Nam.

=> Về mặt kinh tế, mức ROA này cho thấy khả năng sinh lời của tài sản ở mức khá, phản ánh quản trị tài sản hiệu quả và rủi ro thấp. Sự ổn định qua thời gian cho thấy mô hình kinh doanh bền vững, đồng thời là nền tảng quan trọng để cải thiện ROE và duy trì tăng trưởng lợi nhuận trong dài hạn.

Thống kê mô tả cho biến ROE

ds1 %>% summarise(
    mean_roe    = mean(roe, na.rm = TRUE),
    sd_roe      = sd(roe, na.rm = TRUE),
    min_roe     = min(roe, na.rm = TRUE),
    q1_roe      = quantile(roe, 0.25, na.rm = TRUE),
    median_roe  = median(roe, na.rm = TRUE),
    q3_roe      = quantile(roe, 0.75, na.rm = TRUE),
    max_roe     = max(roe, na.rm = TRUE))

Câu lệnh: (1-8) Tạo bảng tóm tắt cho roe về trung bình, độ lệch chuẩn, giá trị nhỏ nhất, phân vị 25%, trung vị, phân vị 75% và giá trị lớn nhất.

Kết quả: ROE trung bình của SHB trong giai đoạn 2014–2023 đạt 0,1129 (≈ 11,29%), nghĩa là mỗi 100 đồng vốn chủ sở hữu tạo ra khoảng 11,29 đồng lợi nhuận sau thuế mỗi năm. Độ lệch chuẩn 0,0366 cho thấy mức biến động vừa phải, phản ánh hiệu quả sinh lời ổn định qua các năm.

Giá trị ROE thấp nhất 0,0690 (6,90%) và cao nhất 0,1801 (18,01%) cho thấy SHB luôn duy trì lợi nhuận dương ngay cả giai đoạn khó khăn, đồng thời có giai đoạn sinh lời vượt trội so với trung bình ngành. Các tứ phân vị lần lượt là Q1 = 0,0822, Median = 0,1066 và Q3 = 0,1384, nghĩa là 50% số năm ROE dao động trong khoảng 8,22%–13,84%, thể hiện phân phối khá cân đối và ổn định.

=> Mức ROE này phản ánh năng lực sử dụng vốn chủ sở hữu hiệu quả, giúp duy trì tăng trưởng lợi nhuận và khả năng tự tài trợ cao. Sự ổn định của ROE qua thời gian cũng cho thấy nền tảng vốn vững chắc, góp phần củng cố niềm tin của cổ đông và nhà đầu tư vào hiệu quả hoạt động của SHB.

Thống kê mô tả cho biến log_total_assets (quy mô tài sản)

ds1 %>% summarise(
    mean_log_ta   = mean(log_total_assets, na.rm = TRUE),
    sd_log_ta     = sd(log_total_assets, na.rm = TRUE),
    min_log_ta    = min(log_total_assets, na.rm = TRUE),
    q1_log_ta     = quantile(log_total_assets, 0.25, na.rm = TRUE),
    median_log_ta = median(log_total_assets, na.rm = TRUE),
    q3_log_ta     = quantile(log_total_assets, 0.75, na.rm = TRUE),
    max_log_ta    = max(log_total_assets, na.rm = TRUE))

Câu lệnh: (1-8) Tạo bảng tóm tắt cho quy mô tài sản về trung bình, độ lệch chuẩn, giá trị nhỏ nhất, phân vị 25%, trung vị, phân vị 75% và giá trị lớn nhất.

Kết quả: Giá trị log trung bình của tổng tài sản đạt 33,46 với độ lệch chuẩn 0,44, cho thấy quy mô tài sản của SHB ổn định và ít biến động qua thời gian. Giá trị nhỏ nhất 32,76 và lớn nhất 34,08 phản ánh mức tăng trưởng đều đặn trong giai đoạn 2014–2023.

Các tứ phân vị Q1 = 33,14, Median = 33,47, Q3 = 33,81 cho thấy phân phối cân đối quanh trung vị, không có biến động đột ngột giữa các năm. Về ý nghĩa kinh tế, điều này chứng tỏ SHB mở rộng quy mô bền vững, tăng trưởng ổn định, không tăng nóng, phản ánh chiến lược phát triển an toàn và quản trị rủi ro hiệu quả trong dài hạn.

Thống kê mô tả cho biến debt_ratio (tỷ lệ nợ)

ds1 %>% summarise(
    mean_debt   = mean(debt_ratio, na.rm = TRUE),
    sd_debt     = sd(debt_ratio, na.rm = TRUE),
    min_debt    = min(debt_ratio, na.rm = TRUE),
    q1_debt     = quantile(debt_ratio, 0.25, na.rm = TRUE),
    median_debt = median(debt_ratio, na.rm = TRUE),
    q3_debt     = quantile(debt_ratio, 0.75, na.rm = TRUE),
    max_debt    = max(debt_ratio, na.rm = TRUE))

Câu lệnh: (1-8) Tạo bảng tóm tắt cho debt_ratio (tỷ lệ nợ) về trung bình, độ lệch chuẩn, giá trị nhỏ nhất, phân vị 25%, trung vị, phân vị 75% và giá trị lớn nhất.

Kết quả: Tỷ lệ nợ trung bình của SHB đạt 0,9373, cho thấy ngân hàng tài trợ khoảng 93,7% tổng tài sản bằng nợ phải trả, tức đòn bẩy tài chính rất cao. Độ lệch chuẩn 0,0142 (≈ 1,42 điểm %) nhỏ phản ánh mức ổn định cao, không có biến động mạnh giữa các năm. Khoảng biến thiên 0,9067–0,9495 cùng các tứ phân vị 0,9319 – 0,9426 – 0,9477 cho thấy phân phối tập trung chặt quanh trung bình, tức cơ cấu vốn ổn định và chính sách sử dụng nợ nhất quán.

=> Tỷ lệ nợ cao giúp tăng hiệu quả đòn bẩy và khả năng mở rộng quy mô tài sản, nhưng cũng làm tăng rủi ro thanh khoản và nghĩa vụ trả nợ nếu môi trường lãi suất bất lợi. Tuy nhiên, độ ổn định của tỷ lệ này cho thấy SHB duy trì được cân bằng hợp lý giữa tăng trưởng và an toàn tài chính.

Thống kê mô tả cho biến cash_ratio (tỷ lệ tiền mặt)

ds1 %>%
  summarise(
    mean_cash   = mean(cash_ratio, na.rm = TRUE),
    sd_cash     = sd(cash_ratio, na.rm = TRUE),
    min_cash    = min(cash_ratio, na.rm = TRUE),
    q1_cash     = quantile(cash_ratio, 0.25, na.rm = TRUE),
    median_cash = median(cash_ratio, na.rm = TRUE),
    q3_cash     = quantile(cash_ratio, 0.75, na.rm = TRUE),
    max_cash    = max(cash_ratio, na.rm = TRUE))

Câu lệnh: (1-8) Tạo bảng tóm tắt cho ash_ratio (tỷ lệ tiền mặt) về trung bình, độ lệch chuẩn, giá trị nhỏ nhất, phân vị 25%, trung vị, phân vị 75% và giá trị lớn nhất.

Kết quả: Tỷ lệ tiền mặt trung bình của SHB đạt 0,0049 (≈ 0,49%), tức chỉ khoảng 0,5% tổng tài sản được nắm giữ dưới dạng tiền mặt và chứng từ có giá trị tương đương tiền. Độ lệch chuẩn 0,0019 cho thấy biến động nhỏ, phản ánh chính sách thanh khoản ổn định qua các năm.

Các tứ phân vị Q1 = 0,0038, Median = 0,0048 và Q3 = 0,0054 cho thấy 50% giá trị nằm trong khoảng hẹp 0,38%–0,54%, tức phân phối tập trung và không có biến động bất thường.

=> Tỷ lệ này phản ánh chiến lược quản trị thanh khoản chặt chẽ, giữ mức tiền mặt tối ưu đủ đáp ứng nhu cầu chi trả ngắn hạn nhưng vẫn tối đa hóa vốn cho hoạt động sinh lời.

Thống kê mô tả cho biến growth_assets (tốc độ tăng trưởng)

ds1 %>%
  summarise(
    mean_growth   = mean(growth_assets, na.rm = TRUE),
    sd_growth     = sd(growth_assets, na.rm = TRUE),
    min_growth    = min(growth_assets, na.rm = TRUE),
    q1_growth     = quantile(growth_assets, 0.25, na.rm = TRUE),
    median_growth = median(growth_assets, na.rm = TRUE),
    q3_growth     = quantile(growth_assets, 0.75, na.rm = TRUE),
    max_growth    = max(growth_assets, na.rm = TRUE))

Câu lệnh: (1-8) Tạo bảng tóm tắt cho growth_assets (tốc độ tăng trưởng tài sản) về trung bình, độ lệch chuẩn, giá trị nhỏ nhất, phân vị 25%, trung vị, phân vị 75% và giá trị lớn nhất.

Kết quả: Tốc độ tăng trưởng tổng tài sản trung bình của SHB đạt 0,1584 (≈ 15,84%), cho thấy ngân hàng mở rộng quy mô ổn định và liên tục qua các năm. Độ lệch chuẩn 0,0494 phản ánh mức biến động vừa phải, tức là tăng trưởng tương đối đều, không xuất hiện đột biến mạnh.

Các tứ phân vị Q1 = 12,99%, Median = 14,29%, Q3 = 21,10% cho thấy 50% giá trị tập trung trong khoảng 13–21%, thể hiện đà tăng trưởng bền vững. Giá trị nhỏ nhất 8,74% và lớn nhất 22,76% khẳng định SHB không có năm nào tăng trưởng âm, duy trì mức mở rộng dương ổn định trong suốt giai đoạn nghiên cứu.

=> Kết quả này phản ánh chiến lược mở rộng tài sản thận trọng nhưng hiệu quả, giúp SHB tăng quy mô hoạt động mà vẫn kiểm soát tốt rủi ro và duy trì tính bền vững trong dài hạn.

Độ lệch (skewness) và độ nhọn (kurtosis)

Câu lệnh dùng để đo đặc điểm phân phối của các biến (ROA, ROE, debt_ratio, cash_ratio): Skewness mức bất đối xứng quanh trung bình và Kurtosis mức tập trung/đỉnh nhọn so với phân phối chuẩn.

library(moments)
vars <- c("roa", "roe", "debt_ratio", "cash_ratio")
res <- dplyr::select(ds1, dplyr::all_of(vars)) |>
  lapply(function(x) c(
    skew = moments::skewness(x, na.rm = TRUE),
    kurt = moments::kurtosis(x, na.rm = TRUE) )) |>
  do.call(what = rbind) |>
  as.data.frame()
round(res, 3)

Câu lệnh: (1) Nạp thư viện tính skewness/kurtosis. (2) Chọn danh sách 4 biến cần đo. (3) Lấy 4 cột đó từ dữ liệu. (4–6) Với mỗi cột: tính skewness và kurtosis (bỏ NA). (7) Ghép kết quả các cột thành bảng nhiều hàng. (8) Đổi sang data.frame. (9) Làm tròn 3 chữ số để hiển thị.

Kết quả: ROA (0,92) và ROE (0,37) lệch phải nhẹ, cho thấy một vài năm SHB đạt lợi nhuận vượt trội. Debt Ratio (–1,14) lệch trái, phản ánh đòn bẩy tài chính cao và ổn định. Cash Ratio (1,13) lệch phải rõ, thể hiện một số năm dự trữ tiền mặt lớn bất thường.

Về độ nhọn, ROA và ROE < 3 cho thấy phân phối phân tán; Debt Ratio ≈ 3 cân bằng; Cash Ratio (4,30) cao nhất, phản ánh tập trung quanh trung bình với vài giá trị cực cao.

=> Nhìn chung, SHB duy trì lợi nhuận tốt trong giai đoạn thuận lợi, sử dụng đòn bẩy ổn định và giữ thanh khoản an toàn, dù còn dư địa tối ưu hóa vốn nhàn rỗi.

So sánh nhóm và giai đoạn

Trung bình ROA theo nhóm quy mô

Mục tiêu là xem xét quy mô tài sản của SHB có ảnh hưởng đến khả năng tạo lợi nhuận trên tổng tài sản hay không, từ đó rút ra xu hướng vận động của hiệu quả kinh doanh theo từng mức độ phát triển.

ds1 %>% group_by(asset_group) %>%
  summarise(mean_roa = mean(roa, na.rm=TRUE), sd_roa = sd(roa, na.rm=TRUE))

Câu lệnh: (1) Nhóm dữ liệu theo asset_group (quy mô tài sản). (2) Tính ROA trung bình (mean_roa) và độ lệch chuẩn (sd_roa) cho từng nhóm, bỏ qua giá trị khuyết.

Nhận xét: Nhóm có quy mô tài sản lớn đạt ROA cao nhất (1,18%), phản ánh khả năng khai thác lợi thế quy mô, tận dụng tốt hơn cơ hội đầu tư và phân bổ nguồn lực. Trong khi đó, nhóm quy mô nhỏ chỉ đạt 0,42%, cho thấy hạn chế về nguồn vốn và hiệu quả vận hành.

Tuy nhiên, độ lệch chuẩn cao hơn ở nhóm lớn (0,0021) cũng cho thấy rủi ro lợi nhuận cao hơn, đặc trưng cho các ngân hàng trong giai đoạn mở rộng nhanh hoặc đa dạng hóa danh mục đầu tư. Điều này hàm ý rằng tăng quy mô giúp tối ưu hóa lợi nhuận, nhưng cần quản lý rủi ro hoạt động và nợ xấu chặt chẽ để duy trì hiệu quả bền vững.

Trung bình ROE theo nhóm quy mô

Phân tích sự thay đổi ROE theo quy mô tài sản nhằm đánh giá xem quy mô hoạt động của SHB có ảnh hưởng tích cực đến khả năng sinh lời vốn chủ sở hữu hay không, từ đó phản ánh hiệu quả sử dụng vốn của cổ đông trong từng giai đoạn phát triển.

ds1 %>% group_by(asset_group) %>% summarise(mean_roe = mean(roe, na.rm=TRUE))

Câu lệnh: (1) Nhóm dữ liệu theo quy mô tài sản. (2) Tính ROE trung bình cho từng nhóm, bỏ qua giá trị khuyết, nhằm so sánh hiệu quả sinh lời trên vốn chủ sở hữu giữa các mức quy mô khác nhau.

Kết quả: Nhóm quy mô tài sản lớn đạt ROE trung bình 15,58%, cao nhất trong ba nhóm, phản ánh khả năng sinh lời vượt trội của SHB khi mở rộng quy mô. Nhóm trung bình đạt 11,16%, thể hiện hiệu quả tài chính ổn định và tích cực, trong khi nhóm nhỏ chỉ đạt 7,17%, cho thấy ở giai đoạn đầu ngân hàng chưa tối ưu hóa việc sử dụng vốn chủ sở hữu.

Xu hướng ROE tăng dần theo quy mô cho thấy lợi thế kinh tế theo quy mô (economies of scale), tức khi SHB mở rộng hoạt động, chi phí bình quân giảm và hiệu quả sử dụng vốn được nâng cao đáng kể, góp phần gia tăng giá trị cho cổ đông.

Trung bình ROA và ROE theo giai đoạn và nhóm đòn bẩy

ds1 %>%
  group_by(period, leverage_group) %>%
  summarise(mean_roa = mean(roa, na.rm = TRUE),
            mean_roe = mean(roe, na.rm = TRUE),
            .groups = "drop")

Câu lệnh: (1) Nhóm dữ liệu theo giai đoạn (period) và mức đòn bẩy (leverage_group). (2) Tính ROA trung bình cho từng nhóm (bỏ NA). (3) Tính ROE trung bình cho từng nhóm (bỏ NA). (4) Bỏ trạng thái nhóm sau khi tính để trả về bảng phẳng

Kết quả: Giai đoạn 2019–2023, SHB đạt ROA và ROE cao hơn rõ rệt so với 2014–2018, cho thấy hiệu quả hoạt động và khả năng sinh lời được cải thiện mạnh. Nhóm đòn bẩy thấp có ROA 1,18% và ROE 15,58%, cao nhất trong các nhóm, chứng tỏ giảm phụ thuộc vào nợ giúp tăng hiệu quả sinh lời. Trong khi đó, nhóm đòn bẩy cao có lợi nhuận thấp hơn do chi phí vốn lớn. Kết quả này khẳng định cơ cấu vốn hợp lý là yếu tố then chốt giúp SHB nâng cao hiệu quả và phát triển bền vững.

So sánh ROA giữa 2 giai đoạn

Việc so sánh này giúp đánh giá sự cải thiện trong hiệu quả sử dụng tài sản, qua đó phản ánh mức độ ổn định và khả năng tạo lợi nhuận của ngân hàng theo thời gian.

ds1 %>% group_by(period) %>% summarise(mean_roa = mean(roa, na.rm=TRUE))

Câu lệnh: (1) Nhóm dữ liệu theo giai đoạn (period), tính ROA trung bình cho từng giai đoạn, bỏ NA và trả về bảng gồm 2 cột: period và mean_roa.

Kết quả: Kết quả cho thấy ROA trung bình của SHB tăng rõ rệt giữa hai giai đoạn. Cụ thể, giai đoạn 2014–2018 chỉ đạt khoảng 0,46%, trong khi 2019–2023 tăng lên gần 0,97%, tức gấp hơn hai lần. Điều này phản ánh sự cải thiện mạnh mẽ về hiệu quả sử dụng tài sản, cho thấy SHB đã vận hành hiệu quả hơn, kiểm soát chi phí tốt hơn và tận dụng quy mô tài sản để gia tăng lợi nhuận.

=> Sự chuyển biến này gắn liền với giai đoạn tăng trưởng mạnh sau tái cơ cấu, khi ngân hàng mở rộng quy mô nhưng vẫn duy trì được hiệu quả sinh lời ổn định.

So sánh ROE giữa các nhóm đòn bẩy

Đánh giá mối quan hệ giữa đòn bẩy tài chính và khả năng sinh lời, qua đó hiểu rõ tác động của việc sử dụng nợ vay đến lợi nhuận của ngân hàng trong giai đoạn 2014–2023.

ds1 %>% group_by(leverage_group) %>%
  summarise(mean_roe = mean(roe, na.rm=TRUE))

Câu lệnh: (1) Nhóm dữ liệu theo leverage_group (mức đòn bẩy). (2) Tính ROE trung bình cho từng nhóm, bỏ NA.

Kết quả: Kết quả cho thấy ROE trung bình thay đổi đáng kể giữa các nhóm đòn bẩy.

Cụ thể, nhóm đòn bẩy thấp có ROE cao nhất (15,58%), phản ánh cấu trúc vốn an toàn giúp tối ưu hóa lợi nhuận trên vốn chủ sở hữu. Nhóm đòn bẩy cao đạt 11,26%, thấp hơn, cho thấy việc sử dụng nợ nhiều làm chi phí tài chính tăng, làm giảm hiệu quả sinh lời. Trong khi đó, nhóm trung bình chỉ đạt 8,09%, thể hiện hiệu quả ở mức thấp hơn kỳ vọng.

=> Kết quả này gợi ý rằng giữ mức đòn bẩy hợp lý giúp SHB cân bằng giữa rủi ro và lợi nhuận, còn đòn bẩy quá cao có thể làm giảm hiệu quả sử dụng vốn trong dài hạn.

So sánh tốc độ tăng trưởng tài sản giữa các nhóm quy mô

Đánh giá xem quy mô tài sản có ảnh hưởng như thế nào đến khả năng mở rộng tổng tài sản của ngân hàng, từ đó rút ra đặc điểm tăng trưởng tương ứng với từng giai đoạn phát triển của SHB.

ds1 %>% group_by(asset_group) %>%
  summarise(mean_growth = mean(growth_assets, na.rm=TRUE))

Câu lệnh: (1) Nhóm dữ liệu theo nhóm quy mô tài sản (asset_group). (2) Tính tốc độ tăng trưởng tài sản trung bình (mean_growth) cho từng nhóm, bỏ qua giá trị NA.

Kết quả: Tốc độ tăng trưởng trung bình của SHB giảm dần khi quy mô nhỏ hơn. Cụ thể, nhóm lớn có mức giảm nhẹ hơn (-0,103), nhóm trung bình giảm sâu hơn (-0,133), và nhóm nhỏ giảm mạnh nhất (-0,160).

Điều này cho thấy các giai đoạn ngân hàng có quy mô nhỏ hơn thường tăng trưởng chậm hơn hoặc chịu biến động cao hơn, trong khi nhóm quy mô lớn duy trì được sự ổn định tốt hơn nhờ lợi thế quy mô và nguồn vốn.

So sánh tỷ lệ nợ debt_ratio giữa các nhóm quy mô

Đánh giá mối quan hệ giữa quy mô hoạt động và mức độ sử dụng đòn bẩy tài chính, từ đó phản ánh chiến lược tài chính và khả năng quản trị rủi ro vốn của ngân hàng.

ds1 %>% group_by(asset_group) %>%
  summarise(mean_debt = mean(debt_ratio, na.rm=TRUE))

Câu lệnh: (1) Phân chia dữ liệu thành ba nhóm theo quy mô tài sản (nhỏ, trung bình, lớn). (2) Tính giá trị trung bình của tỷ lệ nợ (debt_ratio) trong từng nhóm quy mô tài sản.

Kết quả: Tỷ lệ nợ trung bình của SHB dao động từ 0,9190 đến 0,9473, cho thấy mức đòn bẩy tài chính cao và ổn định. Nhóm quy mô tài sản lớn có tỷ lệ nợ thấp nhất (≈91,9%), thể hiện khả năng tự chủ tài chính tốt hơn.Ngược lại, nhóm trung bình và nhỏ có tỷ lệ nợ cao hơn (≈94%), phản ánh sự phụ thuộc lớn hơn vào nguồn vốn vay.

So sánh tỷ lệ nợ (debt_ratio) giữa các nhóm đòn bẩy

So sánh trung bình và độ lệch chuẩn của debt_ratio giữa các nhóm đòn bẩy giúp đánh giá mức nợ điển hình và mức độ biến động của từng nhóm.

ds1 %>% group_by(leverage_group) %>% 
  summarise(mean_debt = mean(debt_ratio, na.rm=TRUE),
            sd_debt = sd(debt_ratio, na.rm=TRUE))

Câu lệnh: (1) Nhóm dữ liệu theo mức đòn bẩy. (2) Tính tỷ lệ nợ trung bình cho từng nhóm (bỏ NA). (3) Tính độ lệch chuẩn của tỷ lệ nợ cho từng nhóm (bỏ NA).

Nhận xét: Nhóm đòn bẩy cao có tỷ lệ nợ trung bình lớn nhất (0,949), thể hiện phụ thuộc mạnh vào nguồn vốn vay. Nhóm trung bình đạt 0,942, thấp hơn nhẹ, còn nhóm thấp chỉ 0,919, phản ánh mức sử dụng nợ thấp nhất. Độ lệch chuẩn nhỏ ở nhóm cao (0,00045) cho thấy mức nợ ổn định, trong khi nhóm thấp biến động lớn hơn (0,0116), chứng tỏ cơ cấu vốn chưa ổn định qua các năm.

=> SHB duy trì đòn bẩy cao ổn định để tối ưu lợi nhuận, nhưng nhóm có nợ thấp lại chịu dao động mạnh hơn về chính sách tài trợ.

Phân tích biến động và rủi ro

Tính hệ số biến thiên

Hệ số biến thiên (CV) đo mức dao động tương đối của dữ liệu so với giá trị trung bình. CV thấp nghĩa là ổn định, còn CV cao cho thấy dao động mạnh và rủi ro lớn hơn.

ds1 %>% summarise(across(c(roa, roe, debt_ratio, cash_ratio),
                   ~sd(.x,na.rm=TRUE)/mean(.x,na.rm=TRUE)))

Câu lệnh: (1) Chọn 4 biến roa, roe, debt_ratio, cash_ratio để tóm tắt. (2) Với mỗi biến, tính CV = sd / mean (đều bỏ NA).

Nhận xét: ROA (0,4892) có mức biến thiên cao nhất, cho thấy hiệu quả sinh lời trên tổng tài sản của SHB dao động đáng kể qua các năm — phản ánh ảnh hưởng của biến động lợi nhuận hoặc thay đổi chiến lược kinh doanh.

ROE (0,3241) cũng dao động tương đối mạnh nhưng ổn định hơn ROA, thể hiện khả năng kiểm soát hiệu quả vốn tốt hơn.

Debt Ratio (0,0151) có CV rất thấp, chứng minh cơ cấu vốn của SHB ổn định và đòn bẩy tài chính được duy trì bền vững.

Cash Ratio (0,3953) biến động khá lớn, cho thấy mức dự trữ tiền mặt thay đổi theo từng giai đoạn để đáp ứng nhu cầu thanh khoản.

Rolling mean (trung bình trượt 3 năm) cho ROA

Tính trung bình trượt 3 năm cho ROA để làm mượt chuỗi thời gian, qua đó nhìn rõ xu hướng dài hạn của hiệu quả hoạt động, thay vì bị nhiễu bởi biến động từng năm (chính sách tín dụng, điều kiện vĩ mô, điều chỉnh vốn).

ds1 <- ds1 %>% arrange(nam) %>%
  mutate(roa_rollmean3 = zoo::rollmean(roa, k = 3, fill = NA, align = "right"))
ds1 %>% dplyr::select(nam, roa, roa_rollmean3)

Câu lệnh: (1) Sắp xếp dữ liệu theo năm. (2) Tạo cột roa_rollmean3 = trung bình trượt 3 năm cho ROA (cửa sổ k=3, canh phải), nên 2014–2015 là NA. (3) Chỉ hiển thị nam, roa, roa_rollmean3.

Kết quả: Hai năm đầu (2014–2015) không có giá trị do chưa đủ 3 quan sát đầu tiên, còn 2023 bị NA vì thiếu dữ liệu năm 2024 — đây là đặc điểm bình thường của trung bình trượt. Giá trị ROA Rolling Mean tăng dần qua thời gian, rõ nhất từ 2018 trở đi, thể hiện xu hướng hiệu quả hoạt động dài hạn cải thiện ổn định. Giai đoạn 2014–2016 ROA còn thấp và biến động nhẹ, trong khi 2019–2023 cho thấy sự tăng trưởng bền vững và hiệu quả hơn.

Phân tích tương quan

Ma trận tương quan

Mục tiêu của ma trận tương quan là đo mức độ và chiều hướng quan hệ tuyến tính giữa các biến tài chính chủ chốt. Kết quả giúp nhận diện các cặp biến liên hệ chặt, phát hiện đa cộng tuyến có thể làm sai lệch ước lượng, và cung cấp căn cứ kinh tế ban đầu về mối quan hệ giữa quy mô, đòn bẩy, thanh khoản, tăng trưởng, hiệu quả sinh lời để định hướng các phân tích/hồi quy tiếp theo

cols <- c("roa","roe","log_total_assets","debt_ratio","cash_ratio","growth_assets")
num_vars <- ds1 |>
  dplyr::select(tidyselect::all_of(cols)) |>
  tidyr::drop_na()
cor(num_vars, use = "pairwise.complete.obs", method = "pearson")

##                         roa        roe log_total_assets debt_ratio
## roa               1,0000000  0,9490924        0,9051046 -0,9358564
## roe               0,9490924  1,0000000        0,9073977 -0,7898412
## log_total_assets  0,9051046  0,9073977        1,0000000 -0,7532667
## debt_ratio       -0,9358564 -0,7898412       -0,7532667  1,0000000
## cash_ratio       -0,7387602 -0,7554188       -0,8904466  0,5915035
## growth_assets    -0,3528139 -0,3801581       -0,3489613  0,3549654
##                  cash_ratio growth_assets
## roa              -0,7387602    -0,3528139
## roe              -0,7554188    -0,3801581
## log_total_assets -0,8904466    -0,3489613
## debt_ratio        0,5915035     0,3549654
## cash_ratio        1,0000000     0,3708455
## growth_assets     0,3708455     1,0000000

Câu lệnh: (1) Xác định danh sách biến dùng để tính tương quan. (2–4) Trích các cột đó từ dữ liệu và loại bỏ giá trị khuyết. (5) Tính ma trận tương quan Pearson (dùng cặp quan sát đầy đủ cho từng cặp biến).

Kết quả:

ROA – ROE: 0,949 → rất cao & dương. Hai thước đo sinh lời nhất quán; khi ROA tăng thì ROE thường tăng cùng chiều.

ROA – log_total_assets: 0,905; ROE – log_total_assets: 0,907 → dương mạnh: quy mô lớn hơn gắn với sinh lời cao hơn (hiệu ứng kinh tế theo quy mô, kỹ năng phân bổ tài sản tốt hơn).

ROA – debt_ratio: −0,936; ROE – debt_ratio: −0,790 → âm mạnh: đòn bẩy cao đi kèm biên lợi nhuận bị bào mòn (chi phí vốn cao, áp lực dự phòng).

ROA – cash_ratio: −0,739; ROE – cash_ratio: −0,755 → âm vừa–mạnh: giữ tiền mặt tương đối nhiều thường kéo giảm hiệu quả sinh lời (tài sản siêu an toàn sinh lợi thấp).

log_total_assets – cash_ratio: −0,890 → rất âm: ngân hàng quy mô lớn có xu hướng nắm giữ tiền mặt tương đối thấp (tận dụng kênh thanh khoản khác hiệu quả hơn).

Growth_assets – ROA/ROE/log_total_assets: −0,353; −0,380; −0,349 (xấp xỉ) → âm nhẹ: năm tăng trưởng tài sản nhanh không nhất thiết đi kèm cải thiện ngay sinh lời hoặc có thể làm loãng biên lợi nhuận ngắn hạn.

Debt_ratio – cash_ratio: ≈ 0,592 → dương vừa: khi đòn bẩy cao, ngân hàng có thể giữ đệm tiền mặt nhỉnh hơn để phòng thủ thanh khoản.

Ma trận hiệp phương sai (Covariance matrix)

cov(num_vars)

##                              roa            roe log_total_assets
## roa               0,000012913206  0,00012347157     0,0012582540
## roe               0,000123471570  0,00131063851     0,0127084224
## log_total_assets  0,001258253983  0,01270842245     0,1496597419
## debt_ratio       -0,000050567190 -0,00042995581    -0,0043817122
## cash_ratio       -0,000005405606 -0,00005568691    -0,0007014298
## growth_assets    -0,000062691028 -0,00068053167    -0,0066753246
##                      debt_ratio      cash_ratio  growth_assets
## roa              -0,00005056719 -0,000005405606 -0,00006269103
## roe              -0,00042995581 -0,000055686907 -0,00068053167
## log_total_assets -0,00438171219 -0,000701429793 -0,00667532460
## debt_ratio        0,00022609197  0,000018110221  0,00026391937
## cash_ratio        0,00001811022  0,000004146176  0,00003733873
## growth_assets     0,00026391937  0,000037338726  0,00244503629

Câu lệnh: (1) Lệnh tính ma trận hiệp phương sai giữa các biến số trong num_vars.

Kết quả: Các hệ số hiệp phương sai dương giữa ROA, ROE và log_total_assets cho thấy ba biến này thường biến động cùng chiều khi quy mô tài sản tăng, hiệu quả sinh lời cũng tăng. Ngược lại, debt_ratio có giá trị âm với ROA và ROE, phản ánh đòn bẩy cao làm giảm hiệu quả sinh lời. Cash_ratio cũng có hiệp phương sai âm với ROA và ROE, cho thấy việc giữ nhiều tiền mặt làm giảm lợi nhuận. Cuối cùng, growth_assets có giá trị dương nhẹ với hầu hết biến sinh lời, biểu thị mối quan hệ mở rộng tài sản gắn với tăng trưởng hiệu quả tài chính.

Trực quan hoá dữ liệu

Xu hướng theo thời gian

ROA xu hướng theo năm

ds1 <- ds1 %>% arrange(nam) %>%
  mutate(lbl = number(roa, accuracy = 0.001, decimal.mark = ","),
    base_y = if_else(row_number() %% 2 == 0,  0.0012, -0.0012),
    dy = case_when(
      nam == 2021 ~  0.0012, 
      nam == 2022 ~  0.0030,  
      nam == 2023 ~ -0.0021,  TRUE        ~  0),
    dx = case_when(
      nam == 2021 ~ -0.12, 
      nam == 2022 ~  0.00,
      nam == 2023 ~  0.12,    
      TRUE        ~  0),
    x_lab = nam + dx,
    y_lab = roa + base_y + dy)
ggplot(ds1, aes(nam, roa)) +
  geom_smooth(method = "loess", se = TRUE, linewidth = 1,
              linetype = "dashed", color = "#1f4e79",
              fill = "#9fb3c8", alpha = 0.25) +
  geom_line(linewidth = 1.4, color = "#2b6cb0") +
  geom_point(size = 3.4, shape = 21, stroke = 1.1, fill = "white", color = "#111111") +
  geom_text(aes(x = x_lab, y = y_lab, label = lbl),
            vjust = 0.5, hjust = 0.5, size = 3.2, family = "Times New Roman", color = "#111111") +
  scale_y_continuous(labels = label_number(accuracy = 0.001, decimal.mark = ","),
                     expand = expansion(mult = c(0.07, 0.18))) +
  scale_x_continuous(breaks = ds1$nam,
                     expand = expansion(mult = c(0.05, 0.08))) +
  labs(title = "ROA CỦA SHB THEO NĂM",
       subtitle = "Đường xanh: ROA thực tế • Đường gạch: xu hướng LOESS (dải xám: sai số chuẩn)",
       x = "Năm", y = "ROA") +
  coord_cartesian(clip = "off") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(
    plot.title    = element_text(face = "bold", size = 12, hjust = 0.5),
    plot.subtitle = element_text(size = 11, color = "#6b7280", hjust = 0.5),
    axis.title    = element_text(size = 12),
    axis.text     = element_text(size = 11),
    panel.grid.minor = element_blank(),
    panel.grid.major.x = element_blank(),
    panel.grid.major.y = element_line(color = "#e5e7eb"),
    plot.margin = margin(10, 24, 12, 12),
    legend.position = "none")

Câu lệnh: (1) Sắp xếp dữ liệu theo năm. (2) Tạo nhãn hiển thị giá trị ROA (làm tròn 3 chữ số). (3) Căn vị trí nhãn xen kẽ để tránh chồng chữ. (4–8) Điều chỉnh vị trí dọc của nhãn (dy) cho các năm đặc biệt. (9–13) Điều chỉnh vị trí ngang (dx) giúp nhãn không đè lên điểm. (14–15) Tạo tọa độ hiển thị nhãn mới (x_lab, y_lab). (16) Khởi tạo biểu đồ ROA theo năm. (17) Thêm đường xu hướng LOESS (gạch xanh nhạt, dải xám = sai số). (18–19) Vẽ đường ROA thực tế và các điểm dữ liệu. (20–21) Gắn nhãn ROA tại từng điểm. (22–24) Định dạng trục và chừa khoảng hiển thị hợp lý. (25–27) Thêm tiêu đề và phụ đề cho biểu đồ. (28–30) Cho phép hiển thị nhãn vượt khung, chọn theme tối giản. (31) Dùng font Times New Roman, cỡ chữ chuẩn. (32–41) Tùy chỉnh thẩm mỹ: căn giữa, màu phụ đề, xóa lưới phụ, thêm lề, ẩn chú giải.

Nhận xét:

ROA của SHB giai đoạn 2014–2023 cho thấy xu hướng tăng rõ rệt, phản ánh hiệu quả sử dụng tài sản được cải thiện theo thời gian. Giai đoạn 2014–2018, ROA duy trì quanh mức thấp (~0,004–0,005), thể hiện giai đoạn tích lũy và mở rộng quy mô. Từ 2019 trở đi, chỉ số tăng nhanh, đạt đỉnh 0,014 năm 2022 trước khi giảm nhẹ còn 0,012 năm 2023.

Đường xu hướng LOESS cho thấy ROA tăng ổn định và bền vững, đặc biệt sau năm 2020 — thời điểm SHB mở rộng hoạt động và nâng cao quản trị hiệu quả. Nhìn chung, mức sinh lời trên tài sản được cải thiện liên tục, cho thấy ngân hàng quản lý tài sản tốt hơn và sử dụng nguồn lực hiệu quả hơn qua từng giai đoạn.

ROE xu hướng theo năm

ggplot(ds1, aes(nam, roe)) +
  geom_ribbon(aes(ymin = 0, ymax = roe), fill = "#CFE4FF", alpha = 0.45) +
  geom_smooth(method = "loess", se = FALSE, color = "#0EA5E9", 
              linewidth = 1, linetype = "dashed") +
  geom_line(linewidth = 1.6, color = "#2563EB") +
  geom_point(shape = 21, size = 3.8, stroke = 1.1, fill = "white", color = "#111111") +
  geom_text_repel(aes(label = scales::number(roe, accuracy = 0.001, decimal.mark = ",")),
                  size = 3.2, family = "Times New Roman", color = "#111111",
                  direction = "y", seed = 123, force = 30, box.padding = 0.3, 
                  point.padding = 0.6, min.segment.length = Inf) +
  scale_y_continuous(labels = scales::label_number(accuracy = 0.001, decimal.mark = ",")) +
  scale_x_continuous(breaks = ds1$nam, expand = expansion(mult = c(0.05, 0.08))) +
  labs(title = "ROE CỦA SHB THEO NĂM",
       subtitle = "Đường xanh: ROE thực tế • Đường gạch: xu hướng LOESS • Vùng xanh nhạt: diện tích tích lũy",
       x = "Năm", y = "ROE") +
  coord_cartesian(clip = "off") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(plot.title = element_text(face = "bold", size = 12, hjust = 0.5),
        plot.subtitle = element_text(size = 11, color = "#6b7280", hjust = 0.5),
        axis.title = element_text(size = 10),
        axis.text = element_text(size = 10),
        panel.grid.minor = element_blank(),
        panel.grid.major.x = element_blank(),
        panel.grid.major.y = element_line(color = "#e5e7eb"),
        legend.position = "none",
        plot.margin = margin(10, 24, 12, 12))

Câu lệnh: (1) Tạo biểu đồ ROE theo năm. (2) Tô vùng tích lũy ROE (màu xanh nhạt). (3) Thêm đường xu hướng LOESS (gạch xanh). (4) Vẽ đường ROE thực tế. (5) Thêm điểm dữ liệu (tròn trắng). (6–9) Gắn nhãn ROE tại từng điểm, căn chỉnh tự động tránh đè nhau. (10–11) Định dạng trục tung – hoành (làm tròn 3 số thập phân, chừa khoảng). (12–15) Thêm tiêu đề, phụ đề, tên trục. (16) Cho phép hiển thị nhãn vượt khung (không bị cắt). (17) Dùng giao diện tối giản, font Times New Roman. (18–26) Tùy chỉnh thẩm mỹ: căn giữa tiêu đề, màu phụ đề, ẩn lưới phụ, chỉ giữ lưới ngang nhạt (#e5e7eb), ẩn chú giải, thêm lề.

Nhận xét

Qua biểu đồ ROE của SHB giai đoạn 2014–2023, có thể thấy tỷ suất lợi nhuận trên vốn chủ sở hữu nhìn chung tăng ổn định, phản ánh hiệu quả sử dụng vốn ngày càng cải thiện.

Giai đoạn 2014–2016, ROE duy trì quanh mức thấp (~0,07), cho thấy hiệu quả sinh lời còn hạn chế trong thời kỳ củng cố hoạt động sau tái cơ cấu.

Từ 2017–2019, ROE tăng dần và đạt 0,131 năm 2019 nhờ mở rộng tín dụng và đa dạng hóa đầu tư.

Giai đoạn 2020–2022, chỉ số tiếp tục tăng mạnh, đạt đỉnh 0,180 năm 2022, thể hiện khả năng sinh lời vượt trội dù chịu tác động dịch bệnh.

Đến 2023, ROE giảm nhẹ xuống 0,146 do chi phí vốn và dự phòng tăng, song vẫn cao hơn đáng kể so với trước 2020.

➡️ Nhìn chung, SHB duy trì xu hướng tăng trưởng ROE bền vững, thể hiện năng lực quản trị và hiệu quả kinh doanh được củng cố, dù cần tiếp tục kiểm soát chi phí và rủi ro tín dụng.

Tỷ lệ nợ/ Tổng tài sản (Debt ratio) theo năm

ggplot(ds1, aes(x = nam, y = 1, fill = debt_ratio)) +
  geom_tile(height = 0.9, width = 0.9) +
  geom_text(aes(label = scales::percent(debt_ratio,
                                        accuracy = 0.1,
                                        decimal.mark = ",")),
            color = "grey20", fontface = "bold", size = 3.6) +
  scale_fill_gradient(low = "#FDE2E4", high = "#F87171",
                      labels = scales::label_percent(accuracy = 0.1,
                                                     decimal.mark = ",")) +
  scale_x_continuous(breaks = ds1$nam) +
  labs(title = "TỶ LỆ NỢ/TỔNG TÀI SẢN — DOT HEATMAP",
       x = "Năm", y = NULL, fill = "Tỷ lệ") +
  theme_minimal(base_family = "Times New Roman") +
  theme(
    axis.text.y  = element_blank(),
    axis.title.y = element_blank(),
    plot.title   = element_text(face = "bold", hjust = 0.5))

Câu lệnh: (1–2) Tạo heatmap thể hiện tỷ lệ nợ theo năm bằng các ô vuông. (3–5) Hiển thị giá trị phần trăm trong từng ô, chữ đậm màu xám đậm, dễ đọc. (6–9) Tô màu chuyển từ hồng nhạt → hồng đậm theo tỷ lệ nợ, định dạng phần trăm có dấu phẩy. (10) Hiển thị năm trên trục X để biểu diễn chuỗi thời gian. (11–12) Thêm tiêu đề, nhãn trục và tên chú giải cho biểu đồ. (13) Dùng giao diện tối giản với font Times New Roman. (14–17) Tùy chỉnh thẩm mỹ: ẩn trục Y, ẩn tiêu đề trục tung, căn giữa và in đậm tiêu đề chính để bố cục gọn, cân đối.

Nhận xét

Giai đoạn 2014–2019: Tỷ lệ dao động nhẹ quanh 94–95%, cho thấy SHB chủ yếu sử dụng nguồn vốn vay và tiền gửi khách hàng làm nguồn tài trợ chính, phản ánh đặc trưng của ngân hàng thương mại với cấu trúc tài sản đòn bẩy cao.

Giai đoạn 2020–2023: Tỷ lệ có xu hướng giảm dần từ 94,9% xuống 90,7% (2022) rồi tăng nhẹ lại 92,1% (2023). Điều này cho thấy SHB giảm phụ thuộc vào nợ, tăng cường vốn chủ sở hữu và cải thiện an toàn tài chính sau giai đoạn tăng trưởng nhanh.

=> Cơ cấu tài chính của SHB vẫn thiên về đòn bẩy cao, song xu hướng gần đây thể hiện chuyển dịch tích cực — giảm dần tỷ lệ nợ, gia tăng nội lực vốn, hướng tới ổn định và bền vững hơn trong dài hạn.

Tỷ lệ tiền mặt (Cash ratio) theo năm

cr <- ds1 %>%
  mutate(
    med = median(cash_ratio, na.rm = TRUE),
    dir = if_else(cash_ratio > med, "Cao hơn trung vị", "Thấp hơn trung vị"),
    lbl = scales::percent(cash_ratio, accuracy = 0.1, decimal.mark = ","),
    y_lab = cash_ratio + if_else(cash_ratio > med, 0.003, -0.003))
cols <- c("Cao hơn trung vị" = "#2563EB", "Thấp hơn trung vị" = "#EF4444")
ggplot(cr, aes(x = nam, y = cash_ratio, color = dir)) +
  geom_hline(aes(yintercept = med), linetype = "longdash", color = "#94A3B8", linewidth = 1) +
  geom_segment(aes(xend = nam, y = med, yend = cash_ratio), linewidth = 1.4) +
  geom_point(size = 4, shape = 21, fill = "white", stroke = 1.1) +
  geom_text(aes(y = y_lab, label = lbl), family = "Times New Roman", size = 3.2) +
  scale_color_manual(values = cols) +
  scale_y_continuous(labels = scales::label_percent(accuracy = 0.1, decimal.mark = ",")) +
  scale_x_continuous(breaks = cr$nam, expand = expansion(mult = c(0.02, 0.02))) +
  labs(
    title = "TỶ LỆ TIỀN MẶT – LOLLIPOP ĐỘ LỆCH SO VỚI TRUNG VỊ",
    subtitle = "Đường gạch: trung vị | Màu xanh: cao hơn trung vị | Màu đỏ: thấp hơn",
    x = "Năm", y = "Cash ratio") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(
    plot.title = element_text(face = "bold", size = 12, hjust = 0.5),
    plot.subtitle = element_text(size = 10, color = "#6b7280", hjust = 0.5),
    legend.position = "top",
    panel.grid.minor = element_blank(),
    panel.grid.major.x = element_blank(),
    panel.grid.major.y = element_line(color = "#E5E7EB"),
    plot.margin = margin(10, 24, 12, 12))

Câu lệnh: (1) Tạo bảng trung gian cr từ ds1. (2–5) Chuẩn bị dữ liệu: tính trung vị, gán nhóm cao/thấp hơn trung vị, tạo nhãn %, chỉnh y_lab để tránh đè. (8) Định nghĩa bảng màu cho 2 nhóm. (9–13) Vẽ lollipop. (14–15) Thang & trục: áp màu, trục y dạng %, trục x theo năm + chừa mép. (16–19) Nhãn biểu đồ: tiêu đề, phụ đề, tên trục x/y. (20) Theme nền: minimal, Times New Roman. (21–23) Tiêu đề & phụ đề: đậm, căn giữa, màu xám cho phụ đề. (24) Chú giải (legend) đặt trên. (25–27) Lưới: tắt lưới phụ & dọc, giữ lưới ngang nhạt để đọc. (28) Canh lề ngoài cho gọn.

Nhận xét

Giai đoạn 2014–2019: Tỷ lệ tiền mặt chủ yếu cao hơn trung vị, thể hiện qua các điểm màu xanh. Đặc biệt năm 2015 vượt trội (≈ 0,9%) cho thấy khả năng thanh khoản ngắn hạn ở mức tốt nhất.

Từ 2020 trở đi: Tỷ lệ chuyển sang thấp hơn trung vị (màu đỏ), giảm dần xuống khoảng 0,2–0,4%, cho thấy xu hướng suy yếu thanh khoản, có thể do tăng chi phí hoặc dòng tiền thuần giảm.

Đường gạch trung vị duy trì ổn định quanh mức 0,45%, giúp dễ nhận thấy các năm lệch nhiều nhất (2015 tăng mạnh, 2023 giảm sâu).

=> Xu hướng tỷ lệ tiền mặt giảm sau 2019, phản ánh mức dự trữ tiền mặt thu hẹp, cần được theo dõi để đảm bảo cân bằng thanh khoản.

Tốc độ tăng trưởng Growth assets theo năm

thr <- 0.10
ga <- ds1 %>%
  filter(!is.na(nam), !is.na(growth_assets)) %>%
  arrange(nam) %>%
  mutate(
    fill_key = if_else(growth_assets >= thr, "ge10", "lt10"),
    lbl = percent(growth_assets, accuracy = 0.1, decimal.mark = ","))
cols_vals <- c(ge10 = "#22C55E", lt10 = "#EF4444")
cols_labs <- c(ge10 = "≥ 10%",   lt10 = "< 10%")
ggplot(ga, aes(nam, growth_assets)) +
  geom_col(aes(fill = fill_key), width = 0.6, color = "white") +
  geom_hline(yintercept = 0,   linewidth = 1,   color = "#94A3B8") +
  geom_hline(yintercept = thr, linewidth = 1.1, color = "#0EA5E9") +
  geom_smooth(method = "loess", se = FALSE, color = "#2563EB",
              linewidth = 1.2, linetype = "dashed") +
  geom_text(aes(label = lbl), vjust = -0.4, size = 3.2,
            family = "Times New Roman", color = "#111111") +
  scale_fill_manual(values = cols_vals, labels = cols_labs,
                    breaks = c("ge10","lt10"), name = NULL, drop = FALSE) +
  scale_y_continuous(labels = label_percent(accuracy = 0.1, decimal.mark = ","),
                     expand = expansion(mult = c(0.08, 0.16))) +
  scale_x_continuous(breaks = ga$nam, expand = expansion(mult = c(0.02, 0.02))) +
  labs(title = "TỐC ĐỘ TĂNG TRƯỞNG TÀI SẢN CỦA SHB THEO NĂM",
       subtitle = "Xanh: ≥ 10% • Đỏ: < 10% • Gạch xanh: LOESS • Đường ngang: 10%",
       x = "Năm", y = "Tốc độ tăng trưởng(%)") +
  coord_cartesian(clip = "off") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(
    plot.title = element_text(face = "bold", size = 14, hjust = 0.5),
    plot.subtitle = element_text(size = 10, color = "#6b7280", hjust = 0.5),
    legend.position = "top")

Câu lệnh: (1) Đặt ngưỡng so sánh tốc độ tăng trưởng 10%. (2–6) Lọc dữ liệu hợp lệ, sắp theo năm, tạo biến nhóm màu (≥10% hoặc <10%) và nhãn phần trăm hiển thị. (7) Xác định màu và nhãn chú giải cho hai nhóm. (8–9) Tạo biểu đồ cột theo năm, tô màu theo nhóm tăng trưởng. (10–13) Kẻ hai đường ngang: mốc 0% và mốc mục tiêu 10% (đậm hơn). (14) Thêm đường xu hướng LOESS (gạch xanh) để mô tả xu thế chung. (15–17) Gắn nhãn phần trăm trên đầu cột, chỉnh vị trí, font, kích cỡ. (18–19) Gán màu và nhãn cho chú giải, ẩn tiêu đề chú giải. (20) Trục tung hiển thị theo định dạng %. (21) Trục hoành hiển thị theo năm, chừa khoảng đệm hai đầu. (22–25) Thêm tiêu đề, phụ đề và tên trục (giải thích ý nghĩa màu, đường gạch, mốc 10%). (26) Cho phép phần tử sát biên không bị cắt. (27–31) Dùng giao diện tối giản, căn giữa tiêu đề & phụ đề, đặt chú giải ở trên.

Nhận xét

Biểu đồ cho thấy tốc độ tăng trưởng tài sản của SHB giai đoạn 2015–2023 dao động nhưng luôn duy trì ở mức dương, phản ánh quy mô tài sản liên tục mở rộng.

Phần lớn các năm, tốc độ tăng trưởng vượt 10%, đặc biệt đạt đỉnh vào 2017 (22,3%) và 2021 (22,8%), cho thấy giai đoạn mở rộng mạnh mẽ. Năm 2022 giảm xuống 8,7%, thấp hơn ngưỡng mục tiêu do ảnh hưởng kinh tế vĩ mô, nhưng 2023 đã phục hồi lên 14,4%, thể hiện xu hướng ổn định trở lại.

=️> Tổng thể, SHB duy trì tăng trưởng tài sản tích cực và bền vững, dù có biến động ngắn hạn, phản ánh năng lực quản trị và mở rộng quy mô hiệu quả.

Phân phối và mô tả

Biểu đồ phân phối ROA

st <- ds1 %>%
  summarise(
    Mean   = mean(roa, na.rm = TRUE),
    Median = median(roa, na.rm = TRUE),
    Min    = min(roa, na.rm = TRUE),
    Max    = max(roa, na.rm = TRUE)) %>%
  mutate(fam = "Times New Roman",
         median_x = 1.06,
         median_y = Median + 0.02)
p <- ggplot(ds1, aes(x = 1, y = roa)) +
  stat_halfeye(adjust = .6, width = .8, .width = 0,
               justification = -.25, point_colour = NA,
               fill = "#93C5FD", color = "#1E3A8A", alpha = .3) +
  geom_boxplot(width = .12, outlier.shape = NA, fill = "white", color = "#1E3A8A") +
  stat_dots(side = "left", justification = 1.1, dotsize = 1.1,
            fill = "#6B7280", color = "#374151")
p +
  geom_point(data = st, inherit.aes = FALSE,
             aes(x = median_x, y = Median),
             size = 3.5, shape = 21, fill = "black", color = "black") +
  geom_segment(data = st, inherit.aes = FALSE,
               aes(x = .75, xend = 1.35, y = Mean, yend = Mean),
               linetype = "longdash", color = "#EF4444", linewidth = .8) +
  annotate("text", x = 1.35, y = st$Mean,
           label = paste0("Mean = ", percent(st$Mean, 0.1, decimal.mark=",")),
           hjust = 0, vjust = -0.6, color = "#EF4444", size = 3.4, family = st$fam) +
  annotate("text", x = st$median_x, y = st$median_y,
           label = paste0("Median = ", percent(st$Median, 0.1, decimal.mark=",")),
           hjust = .5, vjust = 0, color = "#111827", size = 3.4, family = st$fam) +
  annotate("text", x = 1.35, y = st$Max,
           label = paste0("Max = ", percent(st$Max, 0.1, decimal.mark=",")),
           hjust = 0, vjust = -0.6, color = "#2563EB", size = 3.4, family = st$fam) +
  annotate("text", x = 1.35, y = st$Min,
           label = paste0("Min = ", percent(st$Min, 0.1, decimal.mark=",")),
           hjust = 0, vjust = 1.3, color = "#2563EB", size = 3.4, family = st$fam) +
  coord_flip(clip = "off") +
  scale_y_continuous(labels = label_percent(accuracy = 0.1, decimal.mark=",")) +
  labs(
    title    = "PHÂN PHỐI ROA – RAINCLOUD (HALF VIOLIN)",
    subtitle = "Vùng mây: phân phối ROA • Đường đỏ: Mean • Chấm đen: Median • Giá trị cực trị hiển thị bên phải",
    x = NULL, y = "ROA (%)") +
  theme_minimal(base_family = st$fam, base_size = 12) +
  theme(
    plot.title = element_text(face = "bold", hjust = .5, size = 14),
    plot.subtitle = element_text(hjust = .5, size = 10, color = "#4B5563"),
    axis.text.x = element_blank(),
    axis.ticks.x = element_blank(),
    panel.grid.minor = element_blank(),
    panel.grid.major.x = element_blank(),
    plot.margin = margin(14, 80, 12, 20))

Câu lệnh: (1–7) Tóm tắt roa. (8) Tạo ggplot, trục y là ROA. (9–11) Vẽ half-violin biểu thị phân phối. (12–14) Thêm boxplot trắng viền xanh và mây điểm xám ở bên trái. (15) Gọi p để thêm lớp kế tiếp. (16–20) Chấm đen thể hiện Median, thêm đứt đỏ cho Mean. (21–35) Ghi nhãn Mean, Median, Max, Min. (36) Lật trục ngang (coord_flip). (37) Định dạng trục y theo phần trăm (dấu “,”). (38–41) Thêm tiêu đề, phụ đề, tên trục, giải thích ý nghĩa đường và chấm. (42–46) Áp dụng theme minimal, căn giữa tiêu đề, chỉnh cỡ chữ, màu phụ đề. (47–49) Ẩn trục, ticks, lưới nhỏ và lớn. (50) Chỉnh khoảng lề của biểu đồ.

Nhận xét

ROA dao động (0,4–1,4%), trung bình ~0,7%, trung vị ~0,6% ⇒ biên độ vừa phải, lợi nhuận ổn định.

Phân phối lệch phải nhẹ (Mean > Median) do vài năm lợi nhuận cao kéo trung bình lên.

IQR hẹp (~0,55–1,0%) cho thấy hiệu quả hoạt động khá đồng đều. Min ~0,4% cảnh báo biên lợi nhuận thấp, Max ~1,4% thể hiện tiềm năng cải thiện.

→ Hàm ý: cần tối ưu tài sản sinh lợi, kiểm soát rủi ro tín dụng và so sánh hiệu quả với trung vị ngành.

Biểu đồ phân phối ROE

mu  <- mean(ds1$roe, na.rm = TRUE)
med <- median(ds1$roe, na.rm = TRUE)
bw  <- diff(range(ds1$roe, na.rm = TRUE)) / 12
dens <- density(ds1$roe[!is.na(ds1$roe)])
y_top <- max(dens$y, na.rm = TRUE)
ggplot(ds1, aes(x = roe)) +
  geom_histogram(aes(y = after_stat(density)),
                 binwidth = bw, boundary = 0,
                 fill = "#93C5FD", color = "white", alpha = 0.9) +
  geom_density(fill = "#2563EB", color = "#2563EB",
               alpha = 0.18, linewidth = 1.2) +
  geom_vline(xintercept = mu,  color = "#1D4ED8", linewidth = 1.1) +
  geom_vline(xintercept = med, color = "#9CA3AF",
             linewidth = 1, linetype = "longdash") +
  annotate("label", x = med, y = y_top * 4,
           label = paste0("Median: ", percent(med, accuracy = 0.001, decimal.mark = ",")),
           vjust = 0, hjust = 0.5, label.size = 0, color = "#6B7280",
           family = "Times New Roman", size = 6) +  
  scale_x_continuous(labels = label_percent(accuracy = 0.001, decimal.mark = ","),
                     expand = expansion(mult = c(0.02, 0.02))) +
  scale_y_continuous(expand = expansion(mult = c(0.05, 0.30))) +
  labs(
    title = "PHÂN PHỐI ROE",
    subtitle = "Cột: Histogram (mật độ) ● Đường: Kernel density\nVạch: Trung bình (xanh), Trung vị (xám gạch)",
    x = "ROE", y = "Mật độ") +
  coord_cartesian(clip = "off") +
  theme_minimal(base_family = "Times New Roman", base_size = 14) +
  theme(
    plot.title = element_text(face = "bold", size = 18, hjust = 0.5),
    plot.subtitle = element_text(size = 12, color = "#6B7280", hjust = 0.5),
    panel.grid.minor = element_blank(),
    panel.grid.major.x = element_blank(),
    panel.grid.major.y = element_line(color = "#E5E7EB"),
    plot.margin = margin(20, 40, 20, 20))

Câu lệnh: (1–5) Tính Mean, Median, độ rộng bin, mật độ và giá trị đỉnh để định vị nhãn. (6–10) Vẽ histogram và đường mật độ (xanh nhạt, xanh đậm, mờ nhẹ). (11–14) Vẽ hai đường dọc biểu thị Mean (xanh lam) và Median (xám gạch). (15–17) Gắn nhãn “Median: …%” trên đỉnh, font Times New Roman, màu xám. (18–20) Định dạng trục x, y theo phần trăm và mở rộng biên. (21–25) Thêm tiêu đề, phụ đề, nhãn trục, mô tả ý nghĩa. (26–34) Áp dụng theme minimal, phóng to toàn cục, căn giữa, làm gọn lưới, giãn lề để biểu đồ to và rõ.

Nhận xét

Phân phối ROE của SHB tập trung quanh 10–12%, biểu đồ một đỉnh rõ.

Trung vị ~10,66% thấp hơn trung bình, cho thấy lệch phải nhẹ do vài năm ROE cao kéo trung bình lên. Khoảng biến thiên 9–18% phản ánh mức sinh lời ổn định, phân tán vừa phải.

→ ROE hai chữ số, lệch phải nhẹ, hiệu quả ổn định; trung vị ~10–11% phản ánh mức điển hình nên nên ưu tiên xem xét trung vị và độ biến động hơn là giá trị trung bình.

ROA theo nhóm tài sản

st <- ds1 %>% group_by(asset_group) %>%
  summarise(med = median(roa, na.rm=TRUE), .groups="drop")
pal <- setNames(c("#A5B4FC","#93C5FD","#FBCFE8","#FDE68A","#C7F9CC")[seq_along(unique(ds1$asset_group))],
                unique(ds1$asset_group))
ggplot(ds1, aes(x=roa, fill=asset_group)) +
  geom_area(stat="density", alpha=.35, color=NA) +
  geom_vline(data=st, aes(xintercept=med, color=asset_group),
             linetype="longdash", linewidth=.7, show.legend=FALSE) +
  facet_wrap(~asset_group, nrow=1, scales="free_y") +
  scale_fill_manual(values=pal) + scale_color_manual(values=pal) +
  scale_x_continuous(labels=label_percent(.1, decimal.mark=",")) +
  labs(title="PHÂN PHỐI ROA THEO NHÓM QUY MÔ TÀI SẢN",
       x="ROA (%)", y="Mật độ") +
  guides(fill="none") +
  theme_minimal(base_family="Times New Roman", base_size=12) +
  theme(plot.title=element_text(face="bold", hjust=.5),
        panel.grid.minor=element_blank(),
        strip.text=element_text(face="bold"))

Câu lệnh: (1–2) Nhóm theo asset_group, tính trung vị ROA từng nhóm. (3–4) Tạo bảng màu pastel ứng với từng nhóm. (5–6) Vẽ biểu đồ mật độ ROA theo nhóm (vùng màu trong suốt). (7–8) Thêm đường dọc trung vị (gạch đứt). (9) Chia ô theo nhóm, trục y riêng. (10–11) Gán màu và định dạng trục x theo %. (12–13) Thêm tiêu đề và nhãn trục. (14–18) Ẩn chú giải, dùng theme tối giản, font Times New Roman, tiêu đề & nhãn ô in đậm.

Nhận xét:

Nhóm Lớn: Đỉnh ~1,1–1,3%, lệch phải nhẹ, phân phối rộng nhất ⇒ phương sai cao, phản ánh biến động theo chu kỳ tăng trưởng. Median ~1,2% ⇒ hiệu quả vượt trội nhờ lợi thế quy mô và chi phí vốn thấp.

Nhóm Trung bình: Đỉnh ~0,6–0,7%, phân phối hẹp ⇒ ổn định, ít năm cực trị. Median ~0,6% ⇒ hiệu quả trung bình, duy trì ổn định.

Nhóm Nhỏ: Mật độ tập trung ~0,4–0,5%, rất hẹp ⇒ sinh lời thấp, ít biến động nhưng hạn chế về quy mô và chi phí vốn cao.

→ ROA tăng theo quy mô: ~0,5% (Nhỏ), ~0,6% (TB), ~1,2% (Lớn). Nhóm lớn sinh lời cao nhưng biến động hơn; nhóm nhỏ ổn định song kém hiệu quả.

Boxplot ROE theo nhóm đòn bẩy

roe_sum <- ds1 %>% group_by(leverage_group) %>%
  summarise(n    = sum(!is.na(roe)),
    mean = mean(roe, na.rm = TRUE),
    med  = median(roe, na.rm = TRUE),
    sd   = sd(roe, na.rm = TRUE),
    .groups = "drop") %>%
  mutate(se = sd / pmax(sqrt(n), 1),
    ci = qt(0.975, pmax(n - 1, 1)) * se,
    lo = mean - ci,
    hi = mean + ci) %>%
  arrange(mean)
lev_order <- roe_sum$leverage_group
roe_sum$leverage_group <- factor(roe_sum$leverage_group, levels = lev_order)
ggplot(roe_sum, aes(y = leverage_group)) +
  geom_errorbarh(aes(xmin = lo, xmax = hi),
                 height = 0.18, linewidth = 1.2, color = "#3B82F6") +
  geom_point(aes(x = mean), size = 4.2, shape = 21,
             fill = "#0F172A", color = "#0F172A") +
  geom_point(aes(x = med), size = 3.6, shape = 23,
             fill = "#F59E0B", color = "#F59E0B") +
  geom_text(aes(x = hi, label = paste0("Mean ", 
                                       percent(mean, accuracy = 0.1, decimal.mark = ","))),
            nudge_x = 0.02, vjust = 0.5, size = 3.6, family = "Times New Roman", color = "#0F172A") +
  scale_x_continuous(
    labels = label_percent(accuracy = 0.1, decimal.mark = ","),
    expand = expansion(mult = c(0.02, 0.18))) +
  labs(title = "ROE THEO NHÓM ĐÒN BẨY",
    subtitle = "Chấm đen: Mean (kèm 95% CI) | Kim cương cam: Median",
    x = "ROE (%)", y = "Nhóm đòn bẩy") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(
    plot.title  = element_text(face = "bold", hjust = 0.5),
    panel.grid.minor = element_blank())

Câu lệnh: (1–6) Nhóm theo leverage_group, tính số quan sát, Mean, Median, SD; rồi tính sai số chuẩn (SE), khoảng tin cậy 95% (lo, hi) và sắp xếp theo Mean. (7–10) Xác định thứ tự nhóm đòn bẩy theo giá trị trung bình để hiển thị đúng thứ tự trong biểu đồ. (11–22) Vẽ biểu đồ: cột lỗi (CI màu xanh), chấm đen biểu thị Mean, kim cương cam là Median, nhãn chữ “Mean …%” đặt cạnh đầu CI. (23–25) Hiển thị trục x theo %, mở rộng biên nhẹ; thêm tiêu đề, phụ đề và nhãn trục. (30–33) Dùng theme tối giản (Times New Roman), căn giữa tiêu đề, xóa lưới phụ để biểu đồ gọn và rõ.

Nhận xét kết quả

ROE có xu hướng giảm dần theo mức đòn bẩy: Thấp (~15,6%) > Cao (~11,3%) > Trung bình (~8,1%).

Khoảng tin cậy (CI) của từng nhóm không chồng lấn nhiều, cho thấy khác biệt có ý nghĩa.

Nhóm đòn bẩy thấp đạt hiệu quả sinh lời cao nhất nhờ rủi ro tài chính thấp, chi phí lãi vay giảm. Ngược lại, nhóm trung bình và cao có ROE thấp hơn, phản ánh gánh nặng chi phí vốn và hiệu quả sử dụng vốn chưa tương xứng.

→ Khuynh hướng này phù hợp lý thuyết: đòn bẩy cao không luôn giúp tăng ROE, mà chỉ hiệu quả khi doanh nghiệp duy trì lợi nhuận hoạt động đủ lớn để bù chi phí vay

Quan hệ giữa các biến

Mối quan hệ giữa biến log(TA) vs ROA

pts <- ds1 %>%
  transmute(x = as.numeric(log_total_assets),
            y = as.numeric(roa)) %>%
  drop_na()
med_x <- median(pts$x); med_y <- median(pts$y)
ggplot(pts, aes(x, y)) +
  geom_hex(binwidth = c(diff(range(pts$x))/40, diff(range(pts$y))/40)) +
  geom_smooth(method = "loess", se = FALSE, linewidth = 1.5, color = "#2E6FBE") +
  geom_vline(xintercept = med_x, linetype = 2, linewidth = 0.6) +
  geom_hline(yintercept = med_y, linetype = 2, linewidth = 0.6) +
  scale_fill_viridis_c(option = "mako", name = "Số quan sát") +
  scale_y_continuous(labels = label_percent(accuracy = 0.1, decimal.mark = ",")) +
  labs(title = "Quan hệ log(Tổng tài sản) và ROA – Mật độ lục giác (hexbin)",
       x = "log(Tổng tài sản)", y = "ROA (%)") + guides(fill = "none")+
  theme_minimal(base_family = "Times New Roman", base_size = 12)

Câu lệnh: (1–4) Tạo tập dữ liệu gồm log tổng tài sản (trục X) và ROA (trục Y), loại bỏ giá trị thiếu. (5) Tính trung vị cho hai biến để làm mốc so sánh. (6–7) Vẽ biểu đồ mật độ lục giác (hexbin) thể hiện số lượng quan sát theo vùng giá trị. (8) Thêm đường xu hướng hồi quy phi tuyến LOESS để biểu diễn quan hệ tổng thể. (9–10) Vẽ hai đường trung vị (dọc và ngang) để chia không gian thành bốn phần. (11–12) Dùng thang màu viridis (mako) để hiển thị mật độ; trục Y định dạng phần trăm. (13–15) Đặt tiêu đề, nhãn trục, ẩn chú giải, dùng theme tối giản và font Times New Roman.

Nhận xét:

Biểu đồ cho thấy mối quan hệ phi tuyến giữa quy mô tài sản (log) và ROA. Ở nhóm quy mô nhỏ–trung bình, ROA khá ổn định quanh ~0,6%; tuy nhiên, sau ngưỡng log(TS) ≈ 33,5–33,8, ROA tăng mạnh.

Điều này gợi ý sự tồn tại ngưỡng quy mô – khi doanh nghiệp vượt mức tài sản nhất định, hiệu quả sinh lời tăng nhanh, có thể nhờ lợi thế kinh tế theo quy mô hoặc vận hành tối ưu hơn.

Các điểm mật độ tập trung dưới trung vị cho thấy đa số quan sát thuộc nhóm trung bình với ROA thấp, còn các ô phía trên–phải biểu hiện doanh nghiệp lớn có ROA cao hơn.

Khoảng ROA mở rộng ở quy mô lớn phản ánh sự biến động cao (hiện tượng phương sai thay đổi), hàm ý rằng hiệu quả không đồng nhất giữa các nhóm tài sản.

→ Tóm lại: ROA chỉ tăng rõ rệt khi vượt ngưỡng quy mô ~33,5 log(TS); phần lớn quan sát quanh ROA 0,6%, và nhóm lớn có hiệu suất cao nhưng biến động hơn.

Mối quan hệ giữa tỷ lệ nợ (Debt ratio) vs ROA

p_cnt <- ggplot(ds1, aes(debt_ratio, roa)) +
  geom_density_2d_filled(contour_var = "ndensity", alpha = 0.85) +     
  scale_fill_viridis_d(option = "C", direction = -1, name = "Mật độ") +
  geom_point(data = ds1, alpha = 0.6, size = 1.8, inherit.aes = TRUE) +
  geom_smooth(method = "loess", se = FALSE, color = "black") +
  labs(title = "Mối quan hệ giữa tỷ lệ nợ (Debt ratio) vs ROA",
       x = "Tỷ lệ nợ (Debt ratio)", y = "ROA") +
  theme_minimal(base_family = "Times New Roman") +
  theme(plot.title = element_text(face = "bold", hjust = 0.5),
        legend.position = "right")
p_cnt

Câu lệnh: (1–2) Khởi tạo biểu đồ: trục X = debt_ratio, trục Y = ROA. (2) Tô mật độ 2D dạng vùng, độ trong 0.85. (3) Dùng bảng màu viridis option “C”, đảo chiều, đặt tên legend. (4) Rải điểm quan sát chồng lên (alpha 0.6, size 1.8). (5) Vẽ đường LOESS mô tả xu hướng tổng thể, không hiển thị SE, màu đen. (6–7) Đặt tiêu đề và nhãn trục. (8–10) Áp dụng theme minimal (Times New Roman), tiêu đề in đậm, căn giữa, legend bên phải. (11) Hiển thị đối tượng biểu đồ.

Nhận xét:

Biểu đồ cho thấy phần lớn quan sát tập trung ở debt_ratio ≈ 0,94–0,95 với ROA ~0,005–0,006, phản ánh mức đòn bẩy cao nhưng sinh lời thấp và ổn định.

Đường LOESS có độ dốc âm, cho thấy quan hệ nghịch biến – khi tỷ lệ nợ tăng, ROA giảm. Ở vùng debt_ratio thấp hơn (~0,91–0,93), ROA cao hơn (~0,010–0,013), hàm ý nhóm có đòn bẩy thấp hoạt động hiệu quả hơn.

Kết quả phù hợp với lý thuyết chi phí – lợi ích đòn bẩy: nợ giúp tăng lợi nhuận khi thuận lợi nhưng làm giảm ROA trung bình do rủi ro tài chính và chi phí lãi vay cao.

→ Tóm lại, ROA giảm khi đòn bẩy tăng, và chỉ các ngân hàng có tỷ lệ nợ thấp hơn trung bình mới duy trì được hiệu quả sinh lời tốt.

Mối quan hệ giữa ROA và Thu nhập ngoài lãi

pts <- ds1 %>% transmute(x = non_interest_income_share, y = roa) %>% drop_na()
mx <- median(pts$x); my <- median(pts$y)
labp <- label_percent(accuracy = .1, decimal.mark = ",")
tm   <- theme_minimal(base_family = "Times New Roman", base_size = 16)
tb   <- theme(axis.text=element_blank(), axis.ticks=element_blank(),
              panel.grid=element_blank(), plot.margin=margin(0,6,2,6))
vcol <- "#5D5A6D"
p_sc <- ggplot(pts, aes(x,y)) +
  geom_point(color="#7A7D7F", alpha=.45, size=2.2) +
  geom_smooth(method="loess", se=TRUE, color="#3B77D8", fill="#3B77D8",
              linewidth=1.2, alpha=.18) +
  geom_vline(xintercept=mx, linetype=2, linewidth=.7, color=vcol) +
  geom_hline(yintercept=my, linetype=2, linewidth=.7, color=vcol) +
  scale_x_continuous(labels=labp) + scale_y_continuous(labels=labp) +
  labs(x="Tỷ trọng thu nhập ngoài lãi", y="ROA (%)") +
  tm + theme(panel.grid.minor=element_blank(), plot.title.position="plot")
panel_dens <- function(var, fill, v, flip = FALSE){
  g <- ggplot(pts, aes({{var}})) +
    geom_density(fill=fill, alpha=.95) +
    geom_vline(xintercept=v, linetype=2, linewidth=.6, color=vcol) +
    scale_x_continuous(labels=labp) + labs(x=NULL, y=NULL) + tm + tb
  if (flip) g + coord_flip() else g
}
p_top   <- panel_dens(x, "#A6DCEF", mx)
p_right <- panel_dens(y, "#9CC3D5", my, TRUE)
(p_top / (p_sc | p_right) + plot_layout(heights=c(.28,1))) +
  plot_annotation(
    title="THU NHẬP NGOÀI LÃI VÀ ROA",
    theme=theme(plot.title=element_text(family="Times New Roman",
                                        face="bold", size=12, hjust=0),
                plot.margin=margin(6,6,6,6)))

Câu lệnh: (1–3) Tạo tập dữ liệu với hai biến: tỷ trọng thu nhập ngoài lãi (X) và ROA (Y), loại giá trị thiếu. (4) Tính trung vị của X và Y để làm mốc. (5) Định dạng trục theo phần trăm. (6–7) Khai báo theme chung. (8) Chọn một màu trung tính dùng cho các vạch mốc. (9–16) Biểu đồ chính: rải điểm, vẽ đường xu hướng LOESS kèm dải sai số, kẻ hai vạch trung vị (dọc/ngang), đặt nhãn trục %, áp theme và tinh chỉnh lưới/tiêu đề. (17–23) Viết hàm tiện ích tạo biểu đồ mật độ cho một biến (X hoặc Y) kèm vạch trung vị, có tùy chọn xoay cho trục Y; dùng theme đã khai báo. (24–25) Tạo hai biểu đồ biên: phía trên cho X và bên phải cho Y. (26–27) Ghép bố cục: biểu đồ mật độ trên cùng, bên dưới là scatter ở trái và mật độ theo Y ở phải. (28–32) Thêm tiêu đề tổng, canh trái và đặt lề ngoài cho toàn bộ bố cục.

Nhận xét:

Biểu đồ thể hiện mối quan hệ phi tuyến giữa tỷ trọng thu nhập ngoài lãi và ROA.

Ở mức tỷ trọng thấp (khoảng 2–4%), ROA biến động mạnh và thường thấp. Khi tỷ trọng tăng về vùng trung bình 5–6%, ROA cải thiện và tiến gần mức 0%, sau đó gần như đi ngang ở mức cao hơn (7–8%). Điều này cho thấy lợi ích của nguồn thu ngoài lãi chỉ rõ nét khi đạt tỷ trọng vừa phải — vượt quá có thể không làm ROA tăng thêm.

Phân bố mật độ cho thấy phần lớn ngân hàng tập trung ở vùng 5–6%, trong khi ROA dao động quanh 0%, phản ánh hiệu quả sinh lời thấp nhưng không tiêu cực.

Khoảng tin cậy 95% rộng ở vùng tỷ trọng thấp và co lại ở trung tâm, cho thấy độ ổn định cao nhất ở nhóm này.

Tổng thể, kết quả ủng hộ chiến lược đa dạng hóa nguồn thu ở mức vừa phải để cải thiện hiệu quả sinh lời; mở rộng quá mức có thể làm chi phí và rủi ro tăng mà không đem lại lợi ích tương xứng.

Mối quan hệ Cash ratio vs ROA

labn <- scales::label_number(big.mark=".", decimal.mark=",")
ggplot(ds1, aes(cash_ratio, roa)) +
  geom_density_2d_filled(contour_var="ndensity", alpha=.55, show.legend=TRUE) +
  geom_point(size=2.4, alpha=.85, color="#8FC1E3") +
  geom_smooth(method="lm", se=TRUE, linewidth=1.1, color="#E49AB0", fill="#E49AB0", alpha=.15) +
  geom_rug(sides="bl", linewidth=.35, alpha=.7, color="#9ED9CC") +
  ggrepel::geom_text_repel(aes(label=nam), size=3, max.overlaps=18,
                           box.padding=.25, point.padding=.15, seed=1) +
  scale_x_continuous(labels=labn) +
  scale_y_continuous(labels=labn) +
  scale_fill_brewer(name="Mật độ", palette="PuBuGn", direction=1) +
  labs(title="ẢNH HƯỞNG CỦA THANH KHOẢN ĐẾN ROA",
       x="Tỷ lệ tiền mặt (Cash ratio)", y="ROA") +
  theme_minimal(base_family="Times New Roman", base_size=12) +
  theme(plot.title=element_text(face="bold"),
        panel.grid.minor=element_blank(),
        legend.position="right")

Câu lệnh: (1) Tạo hàm định dạng số cho trục theo chuẩn VN. (2–3) Khởi tạo biểu đồ. (4) Rải các điểm quan sát (độ trong cao). (5) Vẽ đường xu hướng tuyến tính (LM) kèm dải sai số. (6) Thêm “rug” ở mép dưới–trái để thấy phân bố dọc theo hai trục. (7–8) Gắn nhãn tên ngân hàng, tự động đẩy nhãn tránh chồng lấn; tinh chỉnh khoảng đệm & hạt giống ngẫu nhiên. (9–10) Định dạng nhãn số cho cả hai trục bằng hàm ở (1). (11) Thang màu cho lớp mật độ, đặt tên chú giải “Mật độ”. (12–13) Đặt tiêu đề và nhãn trục. (14–16) Áp dụng theme tối giản với đúng font; tiêu đề in đậm; ẩn lưới nhỏ. (17) Đặt vị trí chú giải ở bên phải.

Nhận xét:

Biểu đồ cho thấy mối quan hệ nghịch chiều giữa tỷ lệ tiền mặt (Cash ratio) và ROA.

Phần lớn điểm dữ liệu tập trung trong khoảng Cash ratio 0,004–0,005 và ROA 0,007–0,009, thể hiện mức thanh khoản và hiệu quả sinh lời ổn định.

Đường hồi quy có độ dốc âm, cho thấy khi ngân hàng nắm giữ quá nhiều tài sản có tính thanh khoản (tiền mặt, tiền gửi NHNN), hiệu quả sinh lời từ cho vay và đầu tư giảm.

Giai đoạn 2021–2023, các điểm dữ liệu dịch sang ROA cao – Cash ratio thấp, phản ánh việc tối ưu hóa nguồn vốn và giảm tiền nhàn rỗi. Ngược lại, 2015–2018 nằm ở vùng ROA thấp – Cash ratio cao, thể hiện xu hướng giữ thanh khoản phòng thủ, làm giảm lợi nhuận.

Tổng thể, kết quả khẳng định sự đánh đổi giữa an toàn thanh khoản và hiệu quả sinh lời — mức tiền mặt hợp lý giúp cân bằng giữa an toàn tài chính và tối đa hóa lợi nhuận.

Mối quan hệ giữa growth_assets vs ROA

labn <- label_percent(accuracy = 1, decimal.mark = ",")
br   <- scales::pretty_breaks(6)
ggplot(ds1, aes(growth_assets, roa)) +
  stat_density_2d_filled(
    aes(fill = after_stat(level)),
    contour_var = "ndensity", alpha = .85 ) +
  geom_point(size = 1.8, alpha = .55, color = "#3B82F6") +
  stat_ellipse(level = .80, linewidth = .7, linetype = 2, color = "#93C5FD") +
  geom_smooth(method = "lm", se = TRUE, linewidth = .9, color = "#3B82F6", fill = "#93C5FD") +
  geom_text_repel(aes(label = nam), size = 3, max.overlaps = 15, seed = 1) +
  geom_rug(alpha = .25) +
  labs(title = "Tăng trưởng tài sản và ROA", x = "Tăng trưởng tài sản", y = "ROA", fill = "Mật độ") +
  scale_x_continuous(labels = labn, breaks = br) +
  scale_y_continuous(labels = labn, breaks = br) +
  scale_fill_brewer(palette = "YlGnBu", direction = 1, guide = "legend") +
  theme_minimal(base_family = "Times New Roman") +
  theme(plot.title = element_text(face = "bold", size = 12),
        legend.position = "bottom",
        legend.title = element_text(size = 10),
        legend.text  = element_text(size = 9),
        panel.grid.minor = element_blank())

Câu lệnh: (1) Tạo hàm định dạng phần trăm theo chuẩn VN. (2) Tạo hàm sinh mốc đẹp (khoảng 6 mốc) dùng cho cả hai trục. (3) Khởi tạo biểu đồ. (4–6) Nền mật độ 2D (chuẩn hoá theo mật độ), sau đó rải các điểm quan sát. (7) Vẽ ellipse bao phủ ~80% điểm (minh hoạ “đám mây” dữ liệu). (8) Vẽ đường hồi quy tuyến tính. (9–10) Gắn nhãn tên quan sát, tự động đẩy nhãn tránh chồng lấn. (11) Thêm “rug” ở mép trục để thấy phân bố dọc theo hai trục. (12) Đặt tiêu đề, nhãn trục và tên chú giải. (13–14) Định dạng nhãn số cho cả hai trục bằng hàm (1) và mốc (2). (15) Chọn bảng màu cho lớp mật độ, hiển thị chú giải. (16–21) Theme tối giản với font Times New Roman; tiêu đề in đậm; chú giải đặt dưới; chỉnh cỡ chữ tiêu đề/chú giải; ẩn lưới phụ để đồ thị sạch hơn.

Nhận xét:

Biểu đồ cho thấy mối quan hệ đồng biến giữa tăng trưởng tài sản và ROA.

Khi tốc độ tăng tài sản cải thiện (ít âm hơn, gần 0%), ROA tăng rõ, đặc biệt giai đoạn 2021–2022 có hiệu quả cao hơn hẳn so với các năm trước.

Vùng mật độ đậm nhất quanh các năm 2016–2019 phản ánh giai đoạn tăng trưởng chậm, lợi nhuận thấp và ổn định.

Xu hướng này cho thấy khi ngân hàng mở rộng quy mô tài sản hợp lý, hiệu quả sinh lời cũng được cải thiện tương ứng, củng cố vai trò tích cực của tăng trưởng tài sản đối với hiệu quả hoạt động.

Ma trận tương quan

vars <- intersect(c("roa","roe","log_total_assets","debt_ratio","cash_ratio","growth_assets"), names(ds1))
cor_mat <- cor(ds1[, vars], use = "pairwise.complete.obs", method = "pearson")
melted <- reshape2::melt(cor_mat, varnames = c("Var1","Var2"), value.name = "value")
ggplot(melted, aes(Var1, Var2, fill = value)) +
  geom_tile(color = "white", linewidth = 0.4) +
  geom_text(aes(label = scales::number(value, accuracy = 0.01, decimal.mark = ",")), size = 3, color = "black") +
  scale_fill_gradient2(low = "#d73027", mid = "#ffffbf", high = "#4575b4", limits = c(-1,1), midpoint = 0,
                       name = "Hệ số tương quan") +
  labs(title = "Ma trận tương quan giữa các biến tài chính", x = NULL, y = NULL) +
  theme_minimal(base_family = "Times New Roman") +
  theme(plot.title = element_text(face = "bold", size = 12),
        axis.text.x = element_text(angle = 45, hjust = 1),
        panel.grid = element_blank())

Câu lệnh: (1–2) Chọn các biến tài chính, tính ma trận tương quan Pearson giữa chúng, cho phép bỏ qua giá trị thiếu. (3) Chuyển ma trận tương quan sang dạng dài để ggplot đọc được. (4–5) Vẽ ô vuông biểu diễn hệ số tương quan, mỗi ô có viền trắng và hiện giá trị tương quan ngay giữa ô. (6–7) Tạo thang màu: đỏ cho tương quan âm, vàng cho trung tính, xanh cho tương quan dương; giá trị nằm trong khoảng –1 đến 1. (8–9) Đặt tiêu đề, tên trục và chú thích “Hệ số tương quan”. (10–13) Dùng theme tối giản, font Times New Roman, tiêu đề in đậm, chữ trục X nghiêng 45°, bỏ lưới nền phụ để biểu đồ gọn và dễ đọc.

Nhận xét:

1. Mức độ tương quan tổng quát. Ma trận cho thấy tương quan chặt giữa nhiều biến, đặc biệt là ROA, ROE và quy mô tài sản (log_total_assets) có tương quan dương mạnh (xanh đậm), trong khi debt_ratio và cash_ratio có tương quan âm rõ với các biến lợi nhuận (đỏ đậm).

2. Cặp biến lợi nhuận và quy mô. Cặp ROA–ROE (0,95) và ROE–log_total_assets (0,90) cho thấy doanh nghiệp quy mô lớn thường sinh lời cao hơn, phản ánh hiệu quả kinh tế theo quy mô.

3. Cặp biến nợ và khả năng sinh lời. Debt_ratio tương quan âm mạnh với ROA (-0,91) và ROE (-0,74), gợi ý đòn bẩy cao làm giảm lợi nhuận; cash_ratio cũng âm vừa phải (~ -0,70), cho thấy thanh khoản cao đi kèm sinh lời thấp hơn.

4. Tăng trưởng tài sản và các yếu tố khác. Growth_assets tương quan yếu với các biến khác, cho thấy tăng trưởng quy mô không nhất thiết đồng nghĩa với hiệu quả sinh lời hay thay đổi cấu trúc vốn.

So sánh theo nhóm

ROA theo năm, phân tổ theo quy mô

df <- ds1 %>%
  filter(!is.na(growth_assets), !is.na(roa))
bin_n <- 12
df_bins <- df %>%
  group_by(asset_group) %>%
  mutate(bin = cut(growth_assets, breaks = bin_n)) %>%
  group_by(asset_group, bin) %>%
  summarise(x_mean = mean(growth_assets),
            y_mean = mean(roa),
            .groups = "drop")
p <- ggplot(df, aes(growth_assets, roa, color = asset_group)) +
  geom_point(alpha = 0.25, size = 1.7, stroke = 0, show.legend = FALSE) +
  geom_smooth(method = "gam", formula = y ~ s(x, k = 4),
              se = TRUE, linewidth = 1, show.legend = FALSE) +
  geom_point(data = df_bins, aes(x = x_mean, y = y_mean),
             size = 2.4, show.legend = FALSE) +
  geom_line(data = df_bins, aes(x = x_mean, y = y_mean, group = asset_group),
            linewidth = 0.7, show.legend = FALSE) +
  ggrepel::geom_text_repel(
    data = df %>% slice_max(abs(roa), n = 3),
    aes(label = nam),
    size = 3, max.overlaps = 20, show.legend = FALSE) +
  scale_x_continuous(
    labels = label_percent(accuracy = 0.1, decimal.mark = ",", big.mark = "."),
    breaks = breaks_pretty(n = 4),             
    expand = expansion(mult = c(0.02, 0.04))) +
  scale_y_continuous(
    labels = label_percent(accuracy = 0.1, decimal.mark = ",", big.mark = "."),
    breaks = breaks_pretty(n = 5)) +
  scale_color_brewer(palette = "Set2", name = "Nhóm quy mô") +
  labs(
    title = "Tăng trưởng tài sản và ROA theo nhóm quy mô",
    x = "Tốc độ tăng trưởng tài sản (%)",
    y = "ROA (%)") +
  facet_wrap(~ asset_group, nrow = 1, scales = "free_y") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(
    panel.grid.minor = element_blank(),
    legend.position  = "bottom",
    plot.title.position = "plot",
    axis.text.x  = element_text(size = 9),          
    axis.title.x = element_text(margin = margin(t = 8)),
    plot.margin  = margin(6, 10, 14, 10)) +
  guides(x = guide_axis(n.dodge = 1 ))  
p

Câu lệnh: (1–5) Chuẩn bị dữ liệu: lọc thiếu, chọn biến cần, chia theo nhóm quy mô, đặt số “khoảng” tăng trưởng. (6–10) Chia trục tăng trưởng thành các khoảng trong từng nhóm; tính tóm tắt để làm mượt. (11–15) Khởi tạo biểu đồ; rải các điểm gốc để thấy phân bố thực tế theo nhóm. (16–20) Vẽ đường xu hướng mượt và dải tin cậy để nhìn quan hệ tổng thể. (21–25) Đặt các điểm đại diện cho mỗi khoảng và nối lại để tạo “đường mức điển hình” từng nhóm. (26–30) Gắn nhãn vài điểm nổi bật, tự tránh chồng; định dạng trục theo phần trăm, đặt mốc đẹp. (31–35) Tinh chỉnh thang, giãn lề, nhãn trục; đặt tiêu đề/phụ đề ngắn gọn, rõ ý. (36–40) Chia ô biểu đồ theo nhóm quy mô để so sánh song song; kiểm soát chú giải. (41–45) Áp theme tối giản (font, lưới, khoảng cách), căn chỉnh bố cục cuối cùng và vẽ.

Nhận xét:

Biểu đồ cho thấy mối quan hệ không đồng nhất giữa tăng trưởng tài sản và ROA theo quy mô ngân hàng.

Nhóm lớn: ROA giảm khi tăng trưởng tài sản cao hơn, phản ánh chi phí mở rộng hoặc đầu tư lớn làm giảm hiệu quả.

Nhóm nhỏ: Xu hướng giảm mạnh hơn, cho thấy các ngân hàng nhỏ dễ bị ảnh hưởng tiêu cực bởi tăng trưởng nóng.

Nhóm trung bình: Đường dạng chữ U – ROA giảm ở mức tăng trưởng trung bình nhưng cải thiện khi tăng trưởng rất thấp hoặc cao, gợi ý ngưỡng tối ưu.

=️> Kết quả cho thấy tăng trưởng nhanh không luôn đi kèm hiệu quả, đặc biệt với ngân hàng lớn và nhỏ; trong khi nhóm trung bình linh hoạt hơn.

ROE theo giai đoạn

df  <- ds1 %>% filter(!is.na(period), !is.na(roe)) %>% mutate(period = as.factor(period))
s <- df %>% group_by(period) %>%
  summarise(n = n(),
            m = mean(roe, na.rm = TRUE),
            sd = sd(roe, na.rm = TRUE),
            se = sd/sqrt(n),
            t  = qt(0.975, n - 1),
            lo = m - t*se, hi = m + t*se,
            .groups = "drop") %>%
  arrange(m) %>% mutate(period = factor(period, levels = period))
df$period <- factor(df$period, levels = levels(s$period))
mu <- mean(df$roe, na.rm = TRUE)
ggplot(s, aes(x = m, y = period)) +
  geom_point(data = df, aes(x = roe, y = period),
             inherit.aes = FALSE, alpha = .18, size = 1.2,
             position = position_jitter(height = .08), color = "#6C757D") +
  geom_vline(xintercept = mu, linetype = 2, linewidth = .7, color = "#8A8A8A") +
  geom_errorbarh(aes(xmin = lo, xmax = hi), height = 0, linewidth = 3, alpha = .25, color = "#2E6FBE") +
  geom_point(shape = 21, stroke = 1.1, size = 4.5, aes(fill = m), color = "white") +
  scale_fill_gradient(low = "#A6DCEF", high = "#3B77D8", name = "ROE trung bình") +
  scale_x_continuous(labels = label_percent(accuracy = 0.1, decimal.mark = ",")) +
  labs(title = "ROE giữa các giai đoạn", x = "ROE (%)", y = "Giai đoạn") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(panel.grid.minor = element_blank(), legend.position = "right")

Câu lệnh: (1–7): Chuẩn bị & tính toán — lọc dữ liệu, ép period thành factor, tính trung bình, độ lệch chuẩn, SE và khoảng tin cậy 95% cho ROE theo giai đoạn. (8–13): Sắp xếp theo trung bình và đặt lại thứ tự các giai đoạn cho trục y. (16): Tính ROE trung bình toàn kỳ (vẽ đường mốc dọc). (18–21): Vẽ điểm từng năm, vạch trung bình toàn kỳ, thanh sai số và điểm trung bình giai đoạn. (22–23): Tạo thang màu theo ROE trung bình, định dạng trục %. (24–25): Tiêu đề, nhãn, theme tối giản, font Times, ẩn lưới phụ, legend bên phải.

Nhận xét:

Biểu đồ cho thấy ROE tăng rõ rệt từ ~9% (2014–2018) lên ~14–15% (2019–2023), cao hơn trung bình toàn kỳ (~12%). Khoảng tin cậy hai giai đoạn ít chồng lấn, chứng tỏ mức tăng có ý nghĩa thống kê.

Phân tán ROE giai đoạn đầu hẹp và ổn định nhưng thấp, trong khi 2019–2023 rộng và biến động hơn, phản ánh hiệu quả cải thiện song rủi ro cũng tăng.

Tóm lại, ROE tăng nhờ hiệu quả sử dụng vốn và cấu trúc tài sản tốt hơn, nhưng cần duy trì kiểm soát chi phí và rủi ro để ổn định kết quả.

Debt ratio so sánh giữa các nhóm quy mô

mu  <- mean(ds1$debt_ratio, na.rm=TRUE)
dev <- ds1 %>%
  filter(!is.na(debt_ratio), !is.na(asset_group)) %>%
  group_by(asset_group) %>%
  summarise(mean = mean(debt_ratio), .groups="drop") %>%
  mutate(dev = mean - mu,
         label = scales::percent(mean, accuracy=1, decimal.mark=",")) %>%
  arrange(dev) %>%
  mutate(asset_group = factor(asset_group, levels = asset_group))

ggplot(dev, aes(dev, asset_group, fill = dev)) +
  geom_col(width=.6, color="white", linewidth=.4) +
  geom_vline(xintercept=0, color="grey50", linewidth=.7) +
  geom_text(aes(label=label),
            nudge_x=ifelse(dev>=0, .002, -.002),
            hjust=ifelse(dev>=0, 0, 1),
            family="Times New Roman", size=3.5) +
  scale_fill_gradient2(low="#ef4444", mid="#f4f4f5", high="#10b981",
                       midpoint=0, name="Lệch so với trung bình") +
  scale_x_continuous(labels=scales::label_percent(accuracy=1, decimal.mark=",")) +
  labs(title="Tỷ lệ nợ theo nhóm quy mô",
       subtitle="Âm: thấp hơn trung bình | Dương: cao hơn",
       x="Lệch (%)", y=NULL) +
  theme_minimal(base_family="Times New Roman", base_size=15) +
  theme(panel.grid.minor=element_blank(), legend.position="bottom")

Câu lệnh: (1) Tính giá trị trung bình toàn mẫu của tỷ lệ nợ để làm mốc so sánh. (2–7) Lọc dữ liệu, tính trung bình tỷ lệ nợ theo nhóm quy mô, sau đó tạo biến độ lệch so với trung bình và nhãn hiển thị phần trăm. (8–9) Sắp xếp nhóm theo độ lệch và cố định thứ tự factor để giữ đúng thứ tự khi vẽ. (11–12) Vẽ cột thể hiện độ lệch của từng nhóm, có viền trắng. (13) Thêm vạch mốc tại 0 để phân biệt nhóm thấp hơn và cao hơn trung bình. (14) Gắn nhãn phần trăm lên cột, tự căn chỉnh vị trí tùy dấu âm hay dương. (15–18) Áp dụng thang màu phân kỳ (đỏ–xám–xanh) với tâm tại 0, biểu thị mức lệch so với trung bình. (19–20) Định dạng trục X hiển thị theo phần trăm và thêm tiêu đề, phụ đề mô tả ý nghĩa. (21–25) Áp dụng theme tối giản, phông Times New Roman, ẩn lưới phụ và đặt chú giải phía dưới.

Nhận xét:

Biểu đồ cho thấy sự khác biệt về tỷ lệ nợ giữa các nhóm quy mô. Nhóm lớn có tỷ lệ nợ thấp hơn trung bình khoảng 2%, phản ánh mức sử dụng đòn bẩy thận trọng hơn.

Ngược lại, các nhóm nhỏ và trung bình có tỷ lệ nợ cao hơn trung bình, đặc biệt nhóm trung bình vượt nhẹ mức trung bình toàn mẫu.

Điều này gợi ý rằng các ngân hàng quy mô nhỏ và trung bình phụ thuộc nhiều hơn vào vốn vay, trong khi nhóm lớn duy trì cấu trúc vốn an toàn và ổn định hơn.

Growth_assets so sánh giữa các nhóm đòn bẩy

df <- ds1 %>%
  filter(!is.na(leverage_group), !is.na(growth_assets)) %>%
  mutate(leverage_group = factor(leverage_group)) %>%
  group_by(leverage_group) %>%
  mutate(mean_ga = mean(growth_assets)) %>%
  ungroup() %>%
  mutate(leverage_group = fct_reorder(leverage_group, mean_ga))
ggplot(df, aes(growth_assets, leverage_group, fill = after_stat(x))) +
  geom_density_ridges_gradient(
    scale = .9, rel_min_height = .01, color = "white", size = .4,
    quantile_lines = TRUE, quantiles = 2) +
  scale_fill_gradientn(
    colours = c("#FF8C42", "#F6C85F", "#6CCECB", "#4D9DE0", "#A06CD5"),
    name = "Tăng trưởng tài sản (%)") +
  scale_x_continuous(labels = label_percent(accuracy = .1, decimal.mark = ",")) +
  labs(
    title = "TĂNG TRƯỞNG TÀI SẢN THEO NHÓM ĐÒN BẨY",
    x = "Tăng trưởng tài sản (%)", y = "Nhóm đòn bẩy tài chính") +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(panel.grid.minor = element_blank(), legend.position = "right",
        panel.spacing.y = unit(10, "pt"))

Câu lệnh: (1–7) Chuẩn bị dữ liệu: lọc dòng không thiếu, ép leverage_group thành factor, tính trung bình, rồi sắp xếp lại thứ tự theo trung bình đó. (8) Khởi tạo biểu đồ: trục X là tăng trưởng tài sản, trục Y là nhóm đòn bẩy; màu điền lấy từ giá trị X. (9–12) Vẽ đường “ridgeline” mật độ cho từng nhóm: đặt tỉ lệ, chiều cao tối thiểu để loại nhiễu, viền trắng, độ dày nét; kẻ thêm đường phân vị (2 quantiles) để chia dải mật độ. (13–15) Tô màu theo thang liên tục 5 màu tùy biến và đặt tên chú giải “Tăng trưởng tài sản (%)”; định dạng trục X dạng phần trăm với dấu phẩy thập phân. (16–21) Gắn tiêu đề và nhãn trục; dùng theme tối giản với phông Times New Roman, ẩn lưới phụ, đặt chú giải bên phải; tăng khoảng cách dọc giữa các panel để các dải ridgeline thoáng, dễ đọc.

Nhận xét:

Biểu đồ cho thấy mối quan hệ giữa tốc độ tăng trưởng tài sản và mức đòn bẩy tài chính có sự khác biệt rõ rệt giữa các nhóm. Nhóm đòn bẩy thấp có phân bố hẹp, tập trung ở vùng tăng trưởng cao nhất (khoảng –10%), phản ánh khả năng mở rộng tài sản tốt dù sử dụng ít nợ. Trong khi đó, các nhóm đòn bẩy trung bình và cao lại có phân bố nghiêng về vùng tăng trưởng thấp hơn (–15% đến –18%), cho thấy khi tỷ lệ nợ tăng, tốc độ tăng trưởng tài sản có xu hướng chậm lại. Điều này gợi ý rằng đòn bẩy tài chính cao không giúp thúc đẩy tăng trưởng, thậm chí còn có thể kìm hãm hiệu quả mở rộng do chi phí vốn và rủi ro tài chính gia tăng.

Hồ sơ đa biến tài chính chuẩn hoá theo nhóm đòn bẩy

vars    <- intersect(c("roa","roe","log_total_assets","debt_ratio","cash_ratio","growth_assets"),
                     names(ds1))
groupby <- if ("leverage_group" %in% names(ds1)) "leverage_group" else NULL
stopifnot(length(vars) >= 3)
cols <- intersect(c(vars, groupby), names(ds1))
dfz <- as.data.frame(ds1[, cols, drop = FALSE]) %>%
  tidyr::drop_na(all_of(vars)) %>%
  mutate(across(all_of(vars), as.numeric)) %>%
  mutate(across(all_of(vars), ~ (.-mean(., na.rm=TRUE))/sd(., na.rm=TRUE))) %>%
  mutate(id = row_number())   
ord <- vars
if ("roa" %in% vars) {
  cors <- sapply(vars, \(v) cor(dfz[[v]], dfz[["roa"]], use = "pairwise"))
  ord  <- names(sort(cors, decreasing = TRUE))}
df_long <- dfz %>%
  pivot_longer(all_of(vars), names_to = "feature", values_to = "z") %>%
  mutate(feature = factor(feature, levels = ord))
mean_lines <- df_long %>%
  group_by(across(all_of(groupby)), feature) %>%
  summarise(z = mean(z), .groups = "drop")
ggplot(df_long, aes(x = feature, y = z, group = id)) +   
  { if (!is.null(groupby))
      geom_line(aes(color = .data[[groupby]]), alpha = 0.15, linewidth = 0.5)
    else
      geom_line(alpha = 0.15, linewidth = 0.5) } +
  { if (!is.null(groupby))
      geom_line(data = mean_lines,
                aes(x = feature, y = z, group = .data[[groupby]],
                    color = .data[[groupby]]),
                linewidth = 1.6, alpha = 0.95)
    else
      geom_line(data = mean_lines, aes(x = feature, y = z, group = 1),
                linewidth = 1.6, color = "#2E6FBE") } +
  geom_hline(yintercept = 0, linetype = 2, linewidth = 0.6, color = "grey50") +
  scale_y_continuous(breaks = seq(-2, 2, 1)) +
  labs(title = "CẤU TRÚC ĐA BIẾN TÀI CHÍNH",
       subtitle = if (is.null(groupby)) "Đường đậm: trung bình toàn mẫu"
                  else paste(groupby, "· Đường đậm: trung bình từng nhóm"),
       x = "Chỉ tiêu (chuẩn hoá)", y = "Z-score", color = NULL) +
  theme_minimal(base_family = "Times New Roman", base_size = 12) +
  theme(panel.grid.minor = element_blank(),
        axis.text.x = element_text(angle = 20, hjust = 1),
        legend.position = "right")

Câu lệnh: (1–3) Chọn biến tài chính và kiểm tra có cột nhóm hay không. (4–9) Lọc dữ liệu, chuẩn hoá về z-score, thêm id. (10–14) Nếu có roa, sắp xếp biến theo tương quan với nó. (15–17) Chuyển dữ liệu sang dạng dài để vẽ. (18–21) Vẽ đường z-score từng quan sát, tô màu theo nhóm nếu có. (22–33) Vẽ thêm đường trung bình (từng nhóm hoặc toàn mẫu). (34–41) Thêm trục, tiêu đề, chú thích và định dạng theme Times New Roman.

Nhận xét:

Nhóm đòn bẩy thấp có ROA, ROE và quy mô tài sản cao hơn trung bình, nhưng debt_ratio thấp và cash_ratio cũng thấp – phản ánh hiệu quả sinh lời tốt nhưng thanh khoản hạn chế.

Nhóm đòn bẩy cao có debt_ratio và cash_ratio cao, song ROA/ROE thấp hơn, cho thấy họ giữ thanh khoản dự phòng rủi ro nhưng hiệu quả sử dụng vốn kém hơn.

Nhóm trung bình có ROA/ROE thấp nhất, cấu trúc nợ và thanh khoản chỉ ở mức vừa, hàm ý hiệu quả bị kẹt giữa hai cực.

Tăng trưởng tài sản của các nhóm dao động quanh 0, cho thấy khác biệt chủ yếu nằm ở cấu trúc vốn chứ không phải quy mô tăng trưởng.