Trong phân tích hồi quy tuyến tính, việc kiểm tra các giả định của mô hình và xác định các điểm dữ liệu có ảnh hưởng lớn là cực kỳ quan trọng để đảm bảo tính chính xác và độ tin cậy của kết quả. Các công cụ trực quan hóa giúp chúng ta đánh giá nhanh chóng các khía cạnh này. Bài viết này sẽ giới thiệu cách sử dụng Python để tạo ba loại biểu đồ chẩn đoán chính: biểu đồ nhiệt tương quan, biểu đồ kiểm định VIF (Variance Inflation Factor) và bốn biểu đồ phần dư tiêu chuẩn.

Giải thích các Biểu đồ Phần dư Tiêu chuẩn

Bốn biểu đồ phần dư được sử dụng để đánh giá các giả định cơ bản của mô hình hồi quy và phát hiện các vấn đề tiềm ẩn:

1. Biểu đồ Phần dư so với Giá trị Ước tính (Residuals vs Fitted Plot):

   Biểu đồ này hiển thị các giá trị phần dư (sự khác biệt giữa giá trị quan sát và giá trị dự đoán) trên trục tung và các giá trị ước tính (dự đoán của mô hình) trên trục hoành. Mục tiêu là tìm kiếm một sự phân tán ngẫu nhiên của các điểm xung quanh đường 0. Nếu có một mẫu hình rõ ràng (ví dụ: hình nón, đường cong), điều đó có thể cho thấy mối quan hệ phi tuyến tính hoặc phương sai không đồng nhất (heteroscedasticity), vi phạm giả định về tính tuyến tính của mô hình.

2. Biểu đồ Q-Q Chuẩn tắc (Normal Q-Q Plot):

   Biểu đồ Quantile-Quantile chuẩn tắc được dùng để kiểm tra xem các phần dư có tuân theo phân phối chuẩn hay không. Các phần dư được sắp xếp và so sánh với các phân vị lý thuyết của một phân phối chuẩn. Lý tưởng nhất, các điểm dữ liệu nên nằm gần trên một đường thẳng 45 độ. Bất kỳ sự lệch đáng kể nào so với đường thẳng này đều cho thấy các phần dư không được phân phối chuẩn, vi phạm một trong các giả định quan trọng của hồi quy tuyến tính.

3. Biểu đồ Tỷ lệ-Vị trí (Scale-Location Plot / Spread-Location Plot):

   Còn được gọi là biểu đồ Spread-Location, biểu đồ này giúp đánh giá giả định về phương sai đồng nhất (homoscedasticity). Nó vẽ căn bậc hai của giá trị tuyệt đối của phần dư tiêu chuẩn hóa (standardized residuals) trên trục tung so với các giá trị ước tính trên trục hoành. Trong trường hợp lý tưởng, các điểm nên phân tán ngẫu nhiên xung quanh một đường ngang, không có mẫu hình rõ ràng. Sự xuất hiện của một mẫu hình (ví dụ: hình quạt) cho thấy phương sai của phần dư không đồng nhất.

4. Biểu đồ Khoảng cách Cook (Cook’s Distance Plot):

   Biểu đồ Khoảng cách Cook được sử dụng để xác định các điểm dữ liệu có ảnh hưởng lớn (influential observations) đến mô hình hồi quy. Khoảng cách Cook đo lường mức độ thay đổi của các hệ số hồi quy khi một điểm dữ liệu cụ thể được loại bỏ khỏi phân tích. Các điểm có khoảng cách Cook cao hơn ngưỡng nhất định (thường là 4/n hoặc 1, tùy thuộc vào bối cảnh) cần được xem xét kỹ lưỡng vì chúng có thể là ngoại lệ hoặc có ảnh hưởng không cân xứng đến kết quả mô hình.

Triển khai bằng Python

Dưới đây là mã Python được cấu trúc thành các hàm để thực hiện ba loại phân tích này. Chúng ta sẽ sử dụng các thư viện pandas, numpy, matplotlib, seaborn và statsmodels.

import matplotlib.pyplot as plt
import seaborn as sns
import pandas as pd
import numpy as np
import statsmodels.api as sm
import statsmodels.formula.api as smf # Sử dụng cho cú pháp formula

# Cấu hình font cho matplotlib để hỗ trợ tiếng Việt (nếu cần) và hiển thị dấu trừ
plt.rcParams['font.sans-serif'] = ['DejaVu Sans'] # DejaVu Sans thường có sẵn và hỗ trợ Unicode
plt.rcParams['axes.unicode_minus'] = False

def tao_bieu_do_nhiet_tuong_quan(du_lieu_df: pd.DataFrame):
    """
    Tính toán và hiển thị biểu đồ nhiệt hệ số tương quan cho DataFrame.
    """
    print("--- Ma trận tương quan Pearson ---")
    ma_tran_tuong_quan = du_lieu_df.corr().round(2)
    print(ma_tran_tuong_quan)

    plt.figure(figsize=(9, 7))
    sns.heatmap(ma_tran_tuong_quan, annot=True, cmap='viridis', fmt=".2f", linewidths=.5, annot_kws={"size": 10})
    plt.title('Biểu đồ nhiệt hệ số tương quan', fontsize=16)
    plt.xticks(rotation=45, ha='right')
    plt.yticks(rotation=0)
    plt.tight_layout()
    plt.savefig('bieu_do_nhiet_tuong_quan.png')
    plt.show()

def tinh_gia_tri_vif_cho_bien(df_doc_lap: pd.DataFrame, ten_bien: str) -> float:
    """
    Tính toán hệ số VIF (Variance Inflation Factor) cho một biến cụ thể trong tập biến độc lập.
    """
    cac_bien_khac = [col for col in df_doc_lap.columns if col != ten_bien]
    
    if not cac_bien_khac:
        # Nếu không còn biến độc lập nào khác, VIF = 1 (không có đa cộng tuyến)
        return 1.0 
    
    # Tạo công thức hồi quy để dự đoán ten_bien từ các biến độc lập khác
    cong_thuc = f"{ten_bien} ~ {' + '.join(cac_bien_khac)}"
    mo_hinh_hoi_quy = smf.ols(cong_thuc, data=df_doc_lap).fit()
    r_binh_phuong = mo_hinh_hoi_quy.rsquared
    
    # Tránh chia cho 0 nếu R-squared = 1 (đa cộng tuyến hoàn hảo)
    return 1 / (1 - r_binh_phuong) if r_binh_phuong < 1 else float('inf')


def thuc_hien_kiem_dinh_vif(du_lieu_dau_vao: pd.DataFrame, ten_bien_phu_thuoc: str):
    """
    Thực hiện kiểm định VIF cho các biến độc lập trong một mô hình hồi quy 
    và hiển thị kết quả dưới dạng bảng và biểu đồ.
    """
    bien_doc_lap = [col for col in du_lieu_dau_vao.columns if col != ten_bien_phu_thuoc]
    
    if not bien_doc_lap:
        print("Không có biến độc lập nào để tính VIF.")
        return

    # Tạo DataFrame chứa chỉ các biến độc lập
    df_bien_doc_lap = du_lieu_dau_vao[bien_doc_lap]

    ket_qua_vif = {}
    for bien in bien_doc_lap:
        ket_qua_vif[bien] = tinh_gia_tri_vif_cho_bien(df_bien_doc_lap, bien)

    df_vif = pd.DataFrame({'VIF': list(ket_qua_vif.values())}, index=list(ket_qua_vif.keys()))
    df_vif['Dung_sai'] = 1 / df_vif['VIF']
    df_vif = df_vif.sort_values(by='VIF', ascending=False)

    print("\n--- Kết quả kiểm định VIF ---")
    print(df_vif)

    plt.figure(figsize=(10, 6))
    plt.bar(df_vif.index, df_vif['VIF'], color='skyblue')
    plt.axhline(y=5, color='red', linestyle='--', linewidth=2, label="Ngưỡng tham chiếu (VIF = 5)")
    plt.axhline(y=10, color='darkred', linestyle='-.', linewidth=2, label="Ngưỡng tham chiếu (VIF = 10)")
    plt.title("Biểu đồ giá trị VIF của các biến độc lập", fontsize=16)
    plt.xlabel("Biến độc lập", fontsize=12)
    plt.ylabel("Giá trị VIF", fontsize=12)
    plt.xticks(rotation=45, ha='right')
    plt.grid(axis='y', linestyle='--', alpha=0.7)
    plt.legend()
    plt.tight_layout()
    plt.savefig('bieu_do_vif.png')
    plt.show()

def ve_cac_bieu_do_phan_du(ket_qua_mo_hinh: sm.regression.linear_model.OLSResults):
    """
    Vẽ 4 biểu đồ chẩn đoán phần dư chuẩn cho mô hình hồi quy tuyến tính.
    Yêu cầu đầu vào là đối tượng kết quả mô hình từ statsmodels (ví dụ: mo_hinh_ols.fit()).
    """
    gia_tri_du_doan = ket_qua_mo_hinh.fittedvalues
    phan_du_thuc_te = ket_qua_mo_hinh.resid
    phan_du_tieu_chuan = ket_qua_mo_hinh.get_influence().resid_studentized_internal

    fig, axes = plt.subplots(2, 2, figsize=(12, 10))
    fig.suptitle("Biểu đồ chẩn đoán phần dư cho mô hình hồi quy", fontsize=18)

    # 1. Biểu đồ Phần dư so với Giá trị Ước tính (Residuals vs Fitted Plot)
    axes[0, 0].scatter(gia_tri_du_doan, phan_du_thuc_te, alpha=0.7)
    axes[0, 0].set_xlabel('Giá trị ước tính (Fitted Values)', fontsize=12)
    axes[0, 0].set_ylabel('Phần dư (Residuals)', fontsize=12)
    axes[0, 0].set_title('(a) Biểu đồ Phần dư vs. Giá trị ước tính', fontsize=14)
    axes[0, 0].axhline(0, ls='--', color='red', linewidth=1.5)

    # 2. Biểu đồ Q-Q Chuẩn tắc (Normal Q-Q Plot)
    sm.qqplot(phan_du_thuc_te, line='s', ax=axes[0, 1]) # 's' vẽ đường chuẩn tắc (standardized line)
    axes[0, 1].set_title('(b) Biểu đồ Q-Q Chuẩn tắc', fontsize=14)
    axes[0, 1].set_xlabel('Phân vị lý thuyết', fontsize=12)
    axes[0, 1].set_ylabel('Phân vị mẫu', fontsize=12)

    # 3. Biểu đồ Tỷ lệ-Vị trí (Scale-Location Plot / Spread-Location Plot)
    # Trục Y sử dụng căn bậc hai của giá trị tuyệt đối của phần dư tiêu chuẩn hóa
    axes[1, 0].scatter(gia_tri_du_doan, np.sqrt(np.abs(phan_du_tieu_chuan)), alpha=0.7)
    axes[1, 0].set_xlabel('Giá trị ước tính (Fitted Values)', fontsize=12)
    axes[1, 0].set_ylabel('Căn bậc hai của |Phần dư tiêu chuẩn|', fontsize=12)
    axes[1, 0].set_title('(c) Biểu đồ Tỷ lệ-Vị trí', fontsize=14)

    # 4. Biểu đồ Khoảng cách Cook (Cook's Distance Plot)
    influence = ket_qua_mo_hinh.get_influence()
    cooks_dist = influence.cooks_distance[0]
    n_samples = len(cooks_dist)
    
    axes[1, 1].stem(np.arange(n_samples), cooks_dist, basefmt=" ", linefmt="C0-", markerfmt="C0o")
    axes[1, 1].set_xlabel('Chỉ số quan sát', fontsize=12)
    axes[1, 1].set_ylabel('Khoảng cách Cook', fontsize=12)
    axes[1, 1].set_title("(d) Biểu đồ Khoảng cách Cook", fontsize=14)
    
    # Thêm đường tham chiếu cho Khoảng cách Cook (ví dụ: 4/n hoặc 1)
    threshold_4_on_n = 4 / n_samples
    if threshold_4_on_n < 1: # Đảm bảo ngưỡng có ý nghĩa
        axes[1, 1].axhline(y=threshold_4_on_n, color='red', linestyle='--', 
                           label=f'Ngưỡng tham chiếu: 4/n ({threshold_4_on_n:.3f})')
        axes[1, 1].legend(fontsize=10)

    plt.tight_layout(rect=[0, 0.03, 1, 0.95]) # Điều chỉnh layout để tránh tiêu đề chồng chéo
    plt.savefig('bon_bieu_do_phan_du.png')
    plt.show()

Ví dụ minh họa sử dụng các hàm

Để minh họa cách sử dụng các hàm trên, chúng ta sẽ tạo một bộ dữ liệu giả lập và áp dụng các phương pháp phân tích.

# Sử dụng các hàm đã định nghĩa ở trên (hoặc import nếu chúng ở trong một file khác)

# Tạo dữ liệu mẫu
np.random.seed(42)
n_quan_sat = 150
X_gia_lap = np.random.rand(n_quan_sat, 3) * 10 # 3 biến độc lập
df_gia_lap = pd.DataFrame(X_gia_lap, columns=['Chi_phi_quang_cao', 'So_luong_nhan_vien', 'So_gio_lam_viec'])

# Giới thiệu một chút đa cộng tuyến (collinearity) cho ví dụ VIF
df_gia_lap['Chi_phi_phu_tro'] = df_gia_lap['Chi_phi_quang_cao'] * 0.7 + np.random.normal(0, 0.8, n_quan_sat)
df_gia_lap['Muc_do_phuc_tap'] = df_gia_lap['So_luong_nhan_vien'] * 0.5 + df_gia_lap['So_gio_lam_viec'] * 0.3 + np.random.normal(0, 1, n_quan_sat)

# Tạo biến phụ thuộc Y
he_so_that = np.array([1.5, -0.8, 2.0, 0.5, -0.3]) # Hệ số cho các biến (quang_cao, nhan_vien, lam_viec, phu_tro, phuc_tap)
sai_so = np.random.normal(0, 3, n_quan_sat) # Sai số ngẫu nhiên
Y_gia_lap = (
    10 # Hằng số chặn
    + df_gia_lap['Chi_phi_quang_cao'] * he_so_that[0]
    + df_gia_lap['So_luong_nhan_vien'] * he_so_that[1]
    + df_gia_lap['So_gio_lam_viec'] * he_so_that[2]
    + df_gia_lap['Chi_phi_phu_tro'] * he_so_that[3]
    + df_gia_lap['Muc_do_phuc_tap'] * he_so_that[4]
    + sai_so
)

df_du_lieu_day_du = df_gia_lap.copy()
df_du_lieu_day_du['Doanh_thu_Y'] = Y_gia_lap

print("Dữ liệu mẫu (5 hàng đầu):")
print(df_du_lieu_day_du.head())

# 1. Vẽ biểu đồ nhiệt tương quan
print("\n--- Bắt đầu tạo biểu đồ nhiệt tương quan ---")
tao_bieu_do_nhiet_tuong_quan(df_du_lieu_day_du)

# 2. Thực hiện kiểm định VIF
print("\n--- Bắt đầu thực hiện kiểm định VIF ---")
thuc_hien_kiem_dinh_vif(df_du_lieu_day_du, 'Doanh_thu_Y')

# 3. Vẽ 4 biểu đồ phần dư
print("\n--- Bắt đầu tạo 4 biểu đồ phần dư ---")
# Đầu tiên, cần xây dựng một mô hình OLS bằng statsmodels để có đối tượng kết quả mô hình
cong_thuc_mo_hinh = 'Doanh_thu_Y ~ Chi_phi_quang_cao + So_luong_nhan_vien + So_gio_lam_viec + Chi_phi_phu_tro + Muc_do_phuc_tap'
mo_hinh_ols = smf.ols(cong_thuc_mo_hinh, data=df_du_lieu_day_du)
ket_qua_mo_hinh_fit = mo_hinh_ols.fit()

print("\nKết quả tóm tắt mô hình OLS:")
print(ket_qua_mo_hinh_fit.summary())

ve_cac_bieu_do_phan_du(ket_qua_mo_hinh_fit)