Trung bình, trung vị và độ lệch chuẩn trong R

Trung bình và trung vị

Hai đại lượng tóm tắt được dùng nhiều nhất trong thống kê, mỗi cái chỉ cần một lần gọi hàm:

mpg <- mtcars$mpg

print(mean(mpg))
print(median(mpg))

mean() cộng tất cả lại rồi chia cho số lượng. median() sắp xếp các giá trị rồi chọn giá trị ở giữa (hoặc lấy trung bình hai giá trị giữa khi số lượng là chẵn). Với 32 xe trong mtcars, mức tiêu thụ nhiên liệu trung bình khoảng 20,1 mpg và trung vị là 19,2: gần nhau, cho thấy dữ liệu không bị lệch nhiều. Khi chúng chênh nhau nhiều, đó cũng là thông tin; ta sẽ nói đến ở cuối bài.

Một điều khiến ai cũng bất ngờ: nếu vector chứa dù chỉ một NA, cả hai hàm đều trả về NA:

x <- c(4, 8, NA, 6)

print(mean(x))
print(mean(x, na.rm = TRUE))

Đó là chủ ý: R từ chối lặng lẽ giả vờ như giá trị thiếu không tồn tại. na.rm = TRUE nghĩa là "tính trên những giá trị đang có". Gần như mọi hàm tóm tắt trong bài này đều nhận đối số đó. Toàn bộ câu chuyện về cách NA lan truyền nằm trong giá trị thiếu.

Độ lệch chuẩn và phương sai

sd() đo độ phân tán: đại khái, một giá trị điển hình cách trung bình bao xa, theo cùng đơn vị với dữ liệu. var() là bình phương của nó:

mpg <- mtcars$mpg

print(sd(mpg))
print(var(mpg))
print(sd(mpg)^2)

Độ lệch chuẩn khoảng 6 mpg nghĩa là các xe thường nằm trong khoảng chừng 6 mpg quanh mức trung bình 20,1. Vì sd() có cùng đơn vị với dữ liệu, đây là con số bạn báo cáo; var() chủ yếu xuất hiện bên trong các công thức khác.

Đây là chi tiết quan trọng cho bài tập trên lớp: sd() và var() tính thống kê mẫu, tức chia tổng bình phương độ lệch cho n − 1, không phải cho n:

x <- c(2, 4, 4, 4, 5, 5, 7, 9)
n <- length(x)

sample_var <- sum((x - mean(x))^2) / (n - 1)
population_var <- sum((x - mean(x))^2) / n

print(sample_var)      # what var(x) gives you
print(var(x))
print(population_var)  # divide by n instead

Vì sao lại là n − 1? Vì bạn ước lượng trung bình từ chính dữ liệu đó, nên các độ lệch quanh trung bình ước lượng luôn hơi nhỏ hơn thực tế; chia cho n − 1 để bù lại. Vì dữ liệu của bạn hầu như luôn là một mẫu lấy từ một thứ lớn hơn, phiên bản n − 1 là điều bạn cần. Nếu bạn thực sự có toàn bộ tổng thể (mọi học sinh trong lớp, mọi sản phẩm trong danh mục), hãy nhân lại: var(x) * (n - 1) / n.

Sai số chuẩn của trung bình

Độ lệch chuẩn và sai số chuẩn bị nhầm lẫn liên tục, nên hãy tách bạch chúng: sd mô tả dữ liệu, SE mô tả ước lượng trung bình của bạn. R không có hàm se() dựng sẵn, nhưng công thức chỉ một dòng:

mpg <- mtcars$mpg

se <- sd(mpg) / sqrt(length(mpg))
print(se)

Sai số chuẩn giảm khi mẫu lớn lên (thu thập gấp bốn lần dữ liệu thì SE giảm một nửa), vì mẫu lớn hơn xác định trung bình chính xác hơn. sd thì không giảm theo cỡ mẫu; các xe đa dạng đúng như bản chất của chúng, bất kể bạn đo bao nhiêu chiếc. SE là viên gạch nền của khoảng tin cậy, nơi nó thực sự phát huy tác dụng.

summary(): toàn cảnh trong một lần gọi

summary() cho bạn bản tóm tắt năm số cộng thêm trung bình chỉ trong một lần, và nó dùng được trên cả data frame:

print(summary(mtcars$mpg))

print(summary(mtcars[, c("mpg", "wt", "hp")]))

Khi gọi trên một data frame, nó tóm tắt mọi cột: cột số có min, các tứ phân vị, mean và max, còn factor có số đếm cho từng mức. Đây là lệnh đầu tiên nên chạy trên bất kỳ bộ dữ liệu nào vừa tải: những giá trị vô lý (tuổi âm, giá trị lớn nhất là 9999) sẽ lộ ra ngay.

Phân vị, khoảng biến thiên và IQR

quantile() mở rộng trung vị ra mọi điểm cắt:

mpg <- mtcars$mpg

print(quantile(mpg))
print(quantile(mpg, probs = c(0.1, 0.5, 0.9)))

print(range(mpg))
print(IQR(mpg))
print(max(mpg) - min(mpg))

Không truyền đối số thì nó trả về giá trị nhỏ nhất, các tứ phân vị và giá trị lớn nhất. Truyền probs = để lấy các điểm cắt cụ thể: phân vị thứ 10 và thứ 90 ở trên bao quanh vùng chứa phần lớn dữ liệu. IQR() (khoảng cách giữa phân vị thứ 25 và thứ 75) là một thước đo độ phân tán mà, khác với sd(), không bị giá trị ngoại lai kéo lệch. range() trả về min và max thành một cặp.

Mốt: mode() của R KHÔNG làm việc này

Cái bẫy này ai cũng mắc đúng một lần. R có một hàm tên là mode(), và nó chẳng liên quan gì đến thống kê: nó cho biết kiểu lưu trữ của một đối tượng:

x <- c(3, 7, 7, 2, 7, 5, 3)

print(mode(x))   # "numeric" - not what you wanted

counts <- table(x)
print(counts)
print(names(which.max(counts)))

Cách làm cần nhớ: table(x) đếm số lần xuất hiện của mỗi giá trị, which.max() tìm số đếm lớn nhất, và names() lấy ra chính giá trị đó. Lưu ý kết quả là một chuỗi ký tự (tên trong bảng luôn là chuỗi); bọc nó trong as.numeric() nếu bạn cần tính toán. Nếu hai giá trị bằng nhau, which.max() lặng lẽ chỉ trả về giá trị đầu tiên, nên hãy tự xem bảng khi có khả năng hòa.

Trung bình hay trung vị: nên báo cáo cái nào

Trung bình dùng mọi giá trị, đó vừa là điểm mạnh vừa là điểm yếu: chỉ một giá trị cực đoan cũng kéo nó đi. Trung vị chỉ quan tâm đến phần giữa, nên giá trị ngoại lai hầu như không chạm tới nó:

salaries <- c(42000, 45000, 47000, 48000, 51000, 54000, 58000)
print(mean(salaries))
print(median(salaries))

# The CEO joins the dataset:
salaries <- c(salaries, 900000)
print(mean(salaries))
print(median(salaries))

Chỉ thêm một giá trị đã đẩy trung bình từ khoảng 49.300 lên hơn 155.000, một con số không mô tả ai trong dữ liệu cả, trong khi trung vị chỉ nhích từ 48.000 lên 49.500. Đó là lý do thu nhập, giá nhà và thời gian nằm viện được báo cáo bằng trung vị: dữ liệu lệch với đuôi dài khiến trung bình gây hiểu lầm. Với dữ liệu gần đối xứng, hai con số khớp nhau và dùng trung bình là ổn (lại còn hiệu quả hơn về mặt thống kê). Một histogram nhanh sẽ cho bạn biết mình đang ở tình huống nào, và việc so trung bình với trung vị tự nó đã là một phép kiểm tra độ lệch chỉ một dòng.

Những điều cần nhớ

• mean(x) và median(x): thêm na.rm = TRUE khi có giá trị thiếu.
• sd(x) và var(x) tính thống kê mẫu (mẫu số n − 1), đúng thứ bạn cần.
• Sai số chuẩn là sd(x) / sqrt(length(x)): nó đo bạn biết trung bình rõ đến đâu, không phải dữ liệu phân tán ra sao.
• summary() trên một data frame vừa tải là phép kiểm tra nhanh nhất trong R.
• Mốt là names(which.max(table(x))); mode() của R nói về kiểu lưu trữ.
• Dữ liệu lệch hoặc có ngoại lai: báo cáo trung vị. Dữ liệu đối xứng: trung bình là ổn.

Tiếp theo: đo mức độ hai biến cùng biến động, tương quan với cor() và cor.test().