Trong Python, kiểm định t và chi bình phương làm bằng scipy.stats, còn hồi quy tuyến tính và logistic làm bằng statsmodels với cú pháp công thức giống R. Báo cáo đủ thông tin nghĩa là ghi thống kê kiểm định, bậc tự do, p, độ lớn hiệu ứng và khoảng tin cậy, không chỉ ghi p.
Ví dụ trong bài dùng dữ liệu mô phỏng có hạt giống cố định, chạy với scipy 1.18 và statsmodels 0.15, để bạn chạy lại và thấy đúng các con số. Phần cuối nói cách chuyển kết quả sang câu văn theo chuẩn trình bày APA, chuẩn mà nhiều tạp chí khoa học xã hội và tâm lý học yêu cầu.
Kiểm định t: chọn đúng loại và ghi đủ
Welch hay Student
from scipy import stats, rồi r = stats.ttest_ind(a, b, equal_var=False). Tham số equal_var=False chọn kiểm định Welch, không giả định hai nhóm có phương sai bằng nhau; nhiều tài liệu phương pháp khuyên dùng Welch làm mặc định. Mặc định của scipy là equal_var=True, tức kiểm định Student, nên phải khai rõ.
Lấy đủ thông tin
Kết quả có r.statistic, r.pvalue, r.df (với Welch là số lẻ, như 82,85 trong ví dụ), và r.confidence_interval() cho khoảng tin cậy 95% của hiệu hai trung bình. Độ lớn hiệu ứng Cohen d không có sẵn trong hàm này, tự tính bằng hiệu trung bình chia độ lệch chuẩn gộp, hoặc dùng thư viện pingouin.
Mẫu ghép cặp
Đo trước và sau trên cùng người: stats.ttest_rel(truoc, sau). Phi tham số tương ứng: stats.mannwhitneyu() cho hai nhóm độc lập, stats.wilcoxon() cho mẫu ghép cặp.
Có cần kiểm Levene trước rồi mới chọn không
Nhiều giáo trình dạy chạy stats.levene(a, b) trước, phương sai bằng nhau thì dùng Student, khác nhau thì dùng Welch. Cách chọn hai bước này làm sai lệch mức ý nghĩa thực tế, và Welch vẫn cho kết quả gần như Student khi phương sai thật sự bằng nhau. Vì vậy cách gọn hơn là dùng Welch ngay từ đầu và ghi rõ trong bài. Nếu tạp chí hoặc hội đồng quen với kết quả Levene, có thể báo cáo thêm như thông tin mô tả, không dùng làm căn cứ chọn kiểm định.
Chi bình phương và bảng 2×2
Với bảng tần số chéo bang (mảng hoặc kết quả pd.crosstab): c = stats.chi2_contingency(bang). Kết quả có c.statistic, c.pvalue, c.dof và c.expected_freq là bảng tần số kỳ vọng. Kiểm c.expected_freq.min(): khi nhiều ô có tần số kỳ vọng nhỏ, xấp xỉ chi bình phương kém tin cậy; với bảng 2×2, dùng stats.fisher_exact(bang).
Một chi tiết hay bị bỏ qua: với bảng 2×2, chi2_contingency mặc định áp hiệu chỉnh liên tục Yates (correction=True), nên con số có thể khác phần mềm khác hoặc khác phép tính tay. Ghi rõ trong bài có hiệu chỉnh hay không. Độ lớn hiệu ứng cho bảng chéo thường báo bằng Cramér V; scipy có stats.contingency.association(bang, method='cramer').
Đọc bảng summary của hồi quy tuyến tính
import statsmodels.formula.api as smf, rồi m = smf.ols('y ~ x + C(nhom)', data=df).fit() và print(m.summary()). Các phần chính:
| Mục trong bảng | Nghĩa | Có cần báo cáo |
|---|---|---|
coef | Hệ số hồi quy chưa chuẩn hoá | Có |
std err, [0.025 0.975] | Sai số chuẩn và khoảng tin cậy 95% | Có, ít nhất một trong hai |
t và P>|t| | Thống kê t và p của từng hệ số | Có |
R-squared, Adj. R-squared | Tỷ lệ phương sai được giải thích | Có |
F-statistic, Prob (F-statistic) | Kiểm định chung cả mô hình | Thường có |
No. Observations | Số quan sát thật dùng trong mô hình | Có |
Omnibus, Jarque-Bera, Durbin-Watson | Chẩn đoán phần dư | Dùng để kiểm, ít khi báo cáo |
Cond. No. | Số điều kiện, gợi ý đa cộng tuyến hoặc thang đo chênh lệch | Đọc kèm ghi chú cuối bảng |
C(nhom) đánh dấu biến phân loại; chọn nhóm tham chiếu bằng C(nhom, Treatment(reference='B')). Khi lo phương sai sai số thay đổi, dùng fit(cov_type='HC3'): bảng chuyển sang thống kê z và ghi kiểu hiệp phương sai ở phần đầu.
Hồi quy logistic và tỷ số chênh
Biến phụ thuộc phải là số 0 và 1: l = smf.logit('dat ~ x + C(nhom)', data=df).fit(). Hệ số ở thang log odds; tỷ số chênh là np.exp(l.params), khoảng tin cậy của tỷ số chênh là np.exp(l.conf_int()). Bảng summary của logit in Pseudo R-squ., là R bình phương giả kiểu McFadden, không so trực tiếp được với R bình phương của hồi quy tuyến tính.
Nếu fit() báo không hội tụ hoặc cảnh báo tách hoàn toàn (perfect separation), thường có một biến dự báo gần như chia đôi hoàn toàn kết cục. Kiểm bảng chéo của biến đó với kết cục trước khi tìm cách sửa. Số quan sát trong bảng summary của logit cũng cần đối chiếu với cỡ mẫu, như với hồi quy tuyến tính.
Kiểm giả định hồi quy trong statsmodels
Bảng summary không tự kiểm giả định. Phần dư lấy bằng m.resid, giá trị dự báo bằng m.fittedvalues; vẽ phần dư theo giá trị dự báo bằng matplotlib để xem quan hệ có tuyến tính và phương sai có đều không. Kiểm định Breusch–Pagan cho phương sai thay đổi: from statsmodels.stats.diagnostic import het_breuschpagan, rồi het_breuschpagan(m.resid, m.model.exog) trả thống kê LM và p ở hai phần tử đầu.
Đa cộng tuyến: variance_inflation_factor trong statsmodels.stats.outliers_influence tính VIF cho từng cột của ma trận thiết kế m.model.exog. Điểm ảnh hưởng: m.get_influence().cooks_distance[0] cho khoảng cách Cook của từng quan sát. Ghi trong phần phương pháp đã kiểm những gì và xử lý ra sao.
ANOVA một yếu tố và so sánh sau
So trung bình của ba nhóm trở lên: stats.f_oneway(nhom_a, nhom_b, nhom_c) cho F và p nhanh. Muốn bảng ANOVA đầy đủ với tổng bình phương và bậc tự do, chạy mô hình smf.ols('y ~ C(nhom)', data=df).fit() rồi sm.stats.anova_lm(mo_hinh, typ=2), với import statsmodels.api as sm.
Khi F có ý nghĩa và cần biết nhóm nào khác nhóm nào, pairwise_tukeyhsd(df['y'], df['nhom']) trong statsmodels.stats.multicomp in bảng so sánh từng cặp với p đã hiệu chỉnh và khoảng tin cậy. Báo cáo ANOVA theo APA gồm F với hai bậc tự do, p, và độ lớn hiệu ứng như eta bình phương, tính từ bảng anova_lm.
Viết kết quả theo chuẩn APA
Kiểm định t
Một câu đủ ý: nhóm A có điểm trung bình (M = …, SD = …) không khác biệt có ý nghĩa thống kê so với nhóm B (M = …, SD = …), t(82,85) = −1,41, p = 0,161, d = …, khoảng tin cậy 95% của hiệu trung bình [−0,72; 0,12]. Chuẩn APA dùng chữ in nghiêng cho ký hiệu thống kê, và trong bản tiếng Anh bỏ số 0 trước dấu thập phân với những đại lượng không thể vượt quá 1 như p (viết p = .161); bài tiếng Việt thường giữ 0,161 theo quy định trình bày của nơi nộp. Ghi p chính xác thay vì chỉ “p < 0,05”, trừ khi p rất nhỏ thì ghi p < 0,001. Bài tiếng Việt dùng dấu phẩy thập phân hay dấu chấm theo yêu cầu nơi nộp.
Hồi quy
Báo cáo bảng hệ số gồm B, sai số chuẩn hoặc khoảng tin cậy, p cho từng biến; dưới bảng ghi R bình phương, F, bậc tự do, số quan sát và kiểu sai số chuẩn. Với logistic, báo tỷ số chênh kèm khoảng tin cậy.
Đưa số từ Python sang bản thảo
Đừng chép tay từ màn hình. m.summary2().tables[1] trả bảng hệ số dạng DataFrame, ghi ra CSV hay Excel bằng to_csv; câu kết quả có thể dựng bằng f-string với số đã làm tròn, ví dụ f't({r.df:.2f}) = {r.statistic:.2f}, p = {r.pvalue:.3f}'.
Lỗi hay gặp
- Dùng
ttest_indmặc định (Student) mà không biết, rồi ghi trong bài là kiểm định Welch, hoặc ngược lại. - Chỉ ghi “p < 0,05” mà không có thống kê, bậc tự do, độ lớn hiệu ứng.
- Quên rằng
smf.olsbỏ dòng thiếu: số quan sát trong bảng nhỏ hơn cỡ mẫu mô tả ở phần phương pháp mà không giải thích. - Diễn giải tỷ số chênh như tỷ số nguy cơ.
- Chạy hàng chục kiểm định rồi chỉ báo các kiểm định có ý nghĩa, không hiệu chỉnh so sánh bội.
statsmodels.stats.multitest.multipletests()có các phương pháp Bonferroni, Holm, Benjamini–Hochberg.
Câu hỏi thường gặp
statsmodels và scikit-learn khác nhau thế nào khi chạy hồi quy?
statsmodels hướng tới suy diễn thống kê: sai số chuẩn, p, khoảng tin cậy, bảng tóm tắt. scikit-learn hướng tới dự báo và không in các thông tin suy diễn đó. Với bài báo cần kiểm định giả thuyết, dùng statsmodels.
Có cần kiểm phân phối chuẩn trước kiểm định t không?
Với mẫu vừa và lớn, kiểm định t khá bền với lệch chuẩn vừa phải. Xem biểu đồ phân bố quan trọng hơn chạy máy móc kiểm định Shapiro, vốn gần như luôn có ý nghĩa khi mẫu lớn.
Làm sao tính Cohen d trong Python?
Lấy hiệu hai trung bình chia cho độ lệch chuẩn gộp tính từ phương sai và cỡ mẫu hai nhóm. Thư viện pingouin có hàm compute_effsize tính sẵn nhiều loại độ lớn hiệu ứng.
Vì sao kết quả chi bình phương khác SPSS?
Thường do hiệu chỉnh Yates: scipy mặc định áp hiệu chỉnh cho bảng 2×2, còn SPSS in cả hai dòng có và không hiệu chỉnh. Đặt correction=False để so với dòng Pearson của SPSS.
Có cách nào xuất bảng hồi quy đẹp từ statsmodels không?
summary2 trả các bảng dạng DataFrame dễ xuất. Hàm summary_col trong statsmodels.iolib.summary2 ghép nhiều mô hình thành một bảng so sánh, có thể xuất ra dạng chữ, HTML hoặc LaTeX.
