Xóa các ngoại lệ khỏi tính hệ số tương quan

Question

Giả sử chúng ta có hai vectơ số x và y. Hệ số Pearson tương quan giữa x và y được cho bởi

cor (x, y)

Làm thế nào tôi có thể tự động xem xét chỉ là một tập hợp con của x và y trong tính toán (ví dụ 90%) như để tối đa hóa hệ số tương quan?

Answer 1

Nếu bạn thực sự muốn làm điều này (loại bỏ lớn nhất (tuyệt đối) dư), sau đó chúng tôi có thể sử dụng mô hình tuyến tính để ước tính ít nhất giải pháp hình vuông và số dư liên quan và sau đó chọn n% trung bình của dữ liệu. Dưới đây là một ví dụ:

Thứ nhất, tạo ra một số dữ liệu giả:

require(MASS) ## for mvrnorm() 
set.seed(1) 
dat <- mvrnorm(1000, mu = c(4,5), Sigma = matrix(c(1,0.8,1,0.8), ncol = 2)) 
dat <- data.frame(dat) 
names(dat) <- c("X","Y") 
plot(dat) 

Tiếp theo, chúng tôi phù hợp với những mô hình tuyến tính và trích xuất các dư:

res <- resid(mod <- lm(Y ~ X, data = dat)) 

Chức năng quantile() có thể cho chúng ta những yêu cầu số lượng của số dư. Bạn đề nghị giữ lại 90% dữ liệu, vì vậy chúng tôi muốn trên và dưới 0,05 quantiles:

res.qt <- quantile(res, probs = c(0.05,0.95)) 

Chọn những quan sát với dư ở giữa 90% dữ liệu:

want <- which(res >= res.qt[1] & res <= res.qt[2]) 

Chúng ta có thể sau đó hình dung này, với các điểm màu đỏ là những người chúng ta sẽ giữ lại:

plot(dat, type = "n") 
points(dat[-want,], col = "black", pch = 21, bg = "black", cex = 0.8) 
points(dat[want,], col = "red", pch = 21, bg = "red", cex = 0.8) 
abline(mod, col = "blue", lwd = 2) 

Các mối tương quan đối với dữ liệu đầy đủ và tập hợp con được lựa chọn là:

> cor(dat) 
      X   Y 
X 1.0000000 0.8935235 
Y 0.8935235 1.0000000 
> cor(dat[want,]) 
      X   Y 
X 1.0000000 0.9272109 
Y 0.9272109 1.0000000 
> cor(dat[-want,]) 
     X  Y 
X 1.000000 0.739972 
Y 0.739972 1.000000 

Hãy nhận biết rằng ở đây chúng tôi có thể được ném ra dữ liệu hoàn toàn tốt, bởi vì chúng tôi chỉ chọn 5% với dư dương lớn nhất và 5% với lớn nhất tiêu cực. Một cách khác là để chọn 90% với nhỏ tuyệt đối dư:

ares <- abs(res) 
absres.qt <- quantile(ares, prob = c(.9)) 
abswant <- which(ares <= absres.qt) 
## plot - virtually the same, but not quite 
plot(dat, type = "n") 
points(dat[-abswant,], col = "black", pch = 21, bg = "black", cex = 0.8) 
points(dat[abswant,], col = "red", pch = 21, bg = "red", cex = 0.8) 
abline(mod, col = "blue", lwd = 2) 

Với tập hợp con hơi khác nhau này, mối tương quan là hơi thấp:

> cor(dat[abswant,]) 
      X   Y 
X 1.0000000 0.9272032 
Y 0.9272032 1.0000000 

điểm nữa là thậm chí sau đó chúng ta đang ném ra dữ liệu tốt. Bạn có thể muốn nhìn vào khoảng cách Cook như một thước đo sức mạnh của các ngoại lệ, và loại bỏ chỉ những giá trị trên một ngưỡng nhất định của Cook.Wikipedia có thông tin về khoảng cách của Cook và các ngưỡng được đề xuất. Chức năng cooks.distance() có thể được sử dụng để lấy các giá trị từ mod:

> head(cooks.distance(mod)) 
      1   2   3   4   5   6 
7.738789e-04 6.056810e-04 6.375505e-04 4.338566e-04 1.163721e-05 1.740565e-03 

và nếu bạn tính toán ngưỡng (s) gợi ý trên Wikipedia và chỉ loại bỏ những người mà vượt quá ngưỡng. Đối với những dữ liệu này:

> any(cooks.distance(mod) > 1) 
[1] FALSE 
> any(cooks.distance(mod) > (4 * nrow(dat))) 
[1] FALSE 

không ai trong số khoảng cách của Nấu vượt quá ngưỡng đề xuất (. Không gây ngạc nhiên vì cách tôi tạo ra các dữ liệu)

Có nói tất cả điều này, tại sao bạn muốn làm điều này? Nếu bạn đang cố gắng loại bỏ dữ liệu để cải thiện mối tương quan hoặc tạo ra một mối quan hệ đáng kể, điều đó nghe có vẻ hơi giống một chút và giống như dữ liệu nạo vét cho tôi.

Answer 2

Bạn có thể thử bootstrapping dữ liệu của bạn để tìm ra hệ số tương quan cao nhất, ví dụ .:

x <- cars$dist 
y <- cars$speed 
percent <- 0.9   # given in the question above 
n <- 1000    # number of resampling 
boot.cor <- replicate(n, {tmp <- sample(round(length(x)*percent), replace=FALSE); cor(x[tmp], y[tmp])}) 

Và sau khi chạy max(boot.cor). Không được dissapointed nếu tất cả các hệ số tương quan sẽ được tất cả như nhau :)

Answer 3

Điều này có thể đã rõ ràng đối với OP, nhưng chỉ để đảm bảo ... Bạn phải cẩn thận vì cố gắng tối đa hóa mối tương quan có thể thực sự có xu hướng bao gồm ngoại lệ. (@Gavin đã chạm vào điểm này trong câu trả lời/nhận xét của anh ấy). Tôi sẽ là lần đầu tiên xóa các ngoại lệ, rồi tính tương quan. Nói chung, chúng tôi muốn tính toán một mối tương quan mạnh mẽ với các ngoại lệ (và có nhiều phương pháp như vậy trong R).

Chỉ để minh họa điều này một cách đáng kể, chúng ta hãy tạo ra hai vectơ x và y đó là không tương quan:

set.seed(1) 
x <- rnorm(1000) 
y <- rnorm(1000) 
> cor(x,y) 
[1] 0.006401211 

Bây giờ chúng ta hãy thêm một điểm outlier (500,500):

x <- c(x, 500) 
y <- c(y, 500) 

Bây giờ mối tương quan của bất kỳ tập hợp con bao gồm điểm ngoại lệ sẽ gần 100% và sự tương quan của bất kỳ tập hợp con đủ lớn nào không bao gồm ngoại lệ sẽ là gần bằng không. Đặc biệt,

> cor(x,y) 
[1] 0.995741 

Nếu bạn muốn ước tính một "true" tương quan đó không phải là nhạy cảm với giá trị ngoại biên, bạn có thể thử các robust gói:

require(robust) 
> covRob(cbind(x,y), corr = TRUE) 
Call: 
covRob(data = cbind(x, y), corr = TRUE) 

Robust Estimate of Correlation: 
      x   y 
x 1.00000000 -0.02594260 
y -0.02594260 1.00000000 

Bạn có thể chơi xung quanh với các thông số của covRob để quyết định cách cắt dữ liệu. CẬP NHẬT: Ngoài ra còn có rlm (hồi quy tuyến tính mạnh mẽ) trong gói MASS.

Answer 4

Sử dụng method = "spearman" trong cor sẽ mạnh mẽ để gây ô nhiễm và dễ thực hiện vì nó chỉ liên quan đến việc thay thế cor(x, y) bằng cor(x, y, method = "spearman").

Lặp lại phân tích Prasad nhưng sử dụng hệ số tương quan Spearman thay vì chúng ta thấy rằng mối tương quan Spearman là thực sự mạnh mẽ để ô nhiễm ở đây, khôi phục cơ bản không tương quan:

set.seed(1) 

# x and y are uncorrelated 
x <- rnorm(1000) 
y <- rnorm(1000) 
cor(x,y) 
## [1] 0.006401211 

# add contamination -- now cor says they are highly correlated 
x <- c(x, 500) 
y <- c(y, 500) 
cor(x, y) 
## [1] 0.995741 

# but with method = "spearman" contamination is removed & they are shown to be uncorrelated 
cor(x, y, method = "spearman") 
## [1] -0.007270813 

Answer 5

Đây là một khả năng với giá trị ngoại biên bị bắt.Sử dụng một sơ đồ tương tự như Prasad:

library(mvoutlier)  
set.seed(1)  
x <- rnorm(1000)  
y <- rnorm(1000)  
xy <- cbind(x, y)  
outliers <- aq.plot(xy, alpha=0.975) #The documentation/default says alpha=0.025. I think the functions wants 0.975 
cor.plot(x, y)  
color.plot(xy) 
dd.plot(xy) 
uni.plot(xy)  

Trong các câu trả lời khác, 500 bị kẹt ở cuối x và y như một ngoại lệ. Điều đó có thể, hoặc có thể không gây ra vấn đề về bộ nhớ với máy của bạn, vì vậy tôi đã giảm nó xuống còn 4 để tránh điều đó.

x1 <- c(x, 4)  
y1 <- c(y, 4)  
xy1 <- cbind(x1, y1)  
outliers1 <- aq.plot(xy1, alpha=0.975) #The documentation/default says alpha=0.025. I think the functions wants 0.975 
cor.plot(x1, y1)  
color.plot(xy1)  
dd.plot(xy1)  
uni.plot(xy1)  

Dưới đây là những hình ảnh từ x1, y1, xy1 dữ liệu: