加權皮爾森的相關性？

Question

我有一個名爲y的2396x34 double matrix，其中每行（2396）表示由34個連續時間段組成的獨立情況。

我也有一個numeric[34]名爲x，代表了34個連續時間段的單一情況。

目前，我計算每一行之間的相關性在y和x這樣的：

crs[,2] <- cor(t(y),x)

我現在需要的是一個加權相關更換cor功能在上面的語句。權重矢量xy.wt的長度爲34個元素，因此可以爲34個連續時間段中的每一個分配不同的權重。

我發現Weighted Covariance Matrix函數cov.wt，並認爲如果我第一個scale的數據它應該像cor函數一樣工作。實際上你可以指定函數返回一個相關矩陣。不幸的是，我似乎不能以相同的方式使用它，因爲我無法單獨提供我的兩個變量（x和y）。

有沒有人知道我可以用我描述的方式獲得加權相關而不犧牲很多速度？

編輯：也許有些數學函數可以在cor功能之前，爲了得到我正在尋找同樣的結果應用到y。也許如果我乘以xy.wt/sum(xy.wt)每個元素？

編輯＃2我在boot包中發現了另一個函數corr。

corr(d, w = rep(1, nrow(d))/nrow(d)) 

d 
A matrix with two columns corresponding to the two variables whose correlation we wish to calculate. 

w 
A vector of weights to be applied to each pair of observations. The default is equal weights for each pair. Normalization takes place within the function so sum(w) need not equal 1. 

這也不是我所需要的，但它更接近。

編輯＃3 下面是一些代碼來生成的數據我一起工作的類型：

x<-cumsum(rnorm(34)) 
y<- t(sapply(1:2396,function(u) cumsum(rnorm(34)))) 
xy.wt<-1/(34:1) 

crs<-cor(t(y),x) #this works but I want to use xy.wt as weight 

Answer 1

你可以回去了相關的定義。

f <- function(x, y, w = rep(1,length(x))) { 
    stopifnot(length(x) == dim(y)[2]) 
    w <- w/sum(w) 
    # Center x and y, using the weighted means 
    x <- x - sum(x*w) 
    y <- y - apply(t(y) * w, 2, sum) 
    # Compute the variance 
    vx <- sum(w * x * x) 
    vy <- rowSums(w * y * y) # Incorrect: see Heather's remark, in the other answer 
    # Compute the covariance 
    vxy <- colSums(t(y) * x * w) 
    # Compute the correlation 
    vxy/sqrt(vx * vy) 
} 
f(x,y)[1] 
cor(x,y[1,]) # Identical 
f(x, y, xy.wt) 

Answer 2

不幸的是，當y是不止一行的矩陣時，接受的答案是錯誤的。該錯誤是在該行

vy <- rowSums(w * y * y) 

我們希望通過w乘的y列，但是這將通過w的元素，回收需要乘以行。因此

> f(x, y[1, , drop = FALSE], xy.wt) 
[1] 0.103021 

是正確的，因爲在這種情況下，乘法進行逐個元件，這相當於在這裏列方式乘法，但

> f(x, y, xy.wt)[1] 
[1] 0.05463575 

給出了錯誤的答案，由於行嚮明智的倍增。

我們可以校正功能如下

f2 <- function(x, y, w = rep(1,length(x))) { 
    stopifnot(length(x) == dim(y)[2]) 
    w <- w/sum(w) 
    # Center x and y, using the weighted means 
    x <- x - sum(x * w) 
    ty <- t(y - colSums(t(y) * w)) 
    # Compute the variance 
    vx <- sum(w * x * x) 
    vy <- colSums(w * ty * ty) 
    # Compute the covariance 
    vxy <- colSums(ty * x * w) 
    # Compute the correlation 
    vxy/sqrt(vx * vy) 
} 

，並覈對由corr生產的那些從boot包結果：

> res1 <- f2(x, y, xy.wt) 
> res2 <- sapply(1:nrow(y), 
+    function(i, x, y, w) corr(cbind(x, y[i,]), w = w), 
+    x = x, y = y, w = xy.wt) 
> all.equal(res1, res2) 
[1] TRUE 

這本身就給人另一種方式，這個問題可能是解決了。

Answer 3

這裏是計算兩個矩陣之間的加權皮爾遜相關性的概括（而不是矢量和矩陣，如在原來的問題）：

matrix.corr <- function (a, b, w = rep(1, nrow(a))/nrow(a)) 
{ 
    # normalize weights 
    w <- w/sum(w) 

    # center matrices 
    a <- sweep(a, 2, colSums(a * w)) 
    b <- sweep(b, 2, colSums(b * w)) 

    # compute weighted correlation 
    t(w*a) %*% b/sqrt(colSums(w * a**2) %*% t(colSums(w * b**2))) 
} 

使用上面的例子，並從希瑟相關函數，我們可以驗證它：

> sum(matrix.corr(as.matrix(x, nrow=34),t(y),xy.wt) - f2(x,y,xy.wt)) 
[1] 1.537507e-15 

在調用語法而言，這類似於加權cor：

> a <- matrix(c(1,2,3,1,3,2), nrow=3) 
> b <- matrix(c(2,3,1,1,7,3,5,2,8,1,10,12), nrow=3) 
> matrix.corr(a,b) 
    [,1]  [,2] [,3]  [,4] 
[1,] -0.5 0.3273268 0.5 0.9386522 
[2,] 0.5 0.9819805 -0.5 0.7679882 
> cor(a, b) 
    [,1]  [,2] [,3]  [,4] 
[1,] -0.5 0.3273268 0.5 0.9386522 
[2,] 0.5 0.9819805 -0.5 0.7679882