加快成對觀察計數R中

Question

我有一個數據集，其中的每個條目的測量子集隨機丟失：

dat <- matrix(runif(100), nrow=10) 
rownames(dat) <- letters[1:10] 
colnames(dat) <- paste("time", 1:10) 
dat[sample(100, 25)] <- NA 

我很感興趣，在此數據集計算每一行之間的相關性（即AA ，ab，ac，ad，...）。但是，我想通過在結果相關矩陣中將其值設置爲NA來排除少於5個非成對非NA觀測值的相關性。

目前，我這樣做如下：

cor <- cor(t(dat), use = 'pairwise.complete.obs') 
names <- rownames(dat) 
filter <- sapply(names, function(x1) sapply(names, function(x2) 
    sum(!is.na(dat[x1,]) & !is.na(dat[x2,])) < 5)) 
cor[filter] <- NA 

然而，這種操作的實際數據集包含> 1000項非常緩慢。

是否有方法可以基於矢量化方式中的非NA成對觀察數來過濾單元格，而不是在嵌套循環中？

Answer 1

您可以使用矩陣方法計算非NA成對觀測值的數量。

讓我們使用這個數據生成代碼。我使數據更大並增加了更多的NAs。

nr = 1000; 
nc = 900; 
dat = matrix(runif(nr*nc), nrow=nr) 
rownames(dat) = paste(1:nr) 
colnames(dat) = paste("time", 1:nc) 
dat[sample(nr*nc, nr*nc*0.9)] = NA 

然後你過濾代碼正在85秒

tic = proc.time() 
names = rownames(dat) 
filter = sapply(names, function(x1) sapply(names, function(x2) 
    sum(!is.na(dat[x1,]) & !is.na(dat[x2,])) < 5)); 
toc = proc.time(); 
show(toc-tic); 
# 85.50 seconds 

我的版本創建與原始數據值1用於非NAS的矩陣。然後使用矩陣乘法計算配對非NAs的數量。它跑了幾分之一秒。

tic = proc.time() 
NAmat = matrix(0, nrow = nr, ncol = nc) 
NAmat[ !is.na(dat) ] = 1; 
filter2 = (tcrossprod(NAmat) < 5) 
toc = proc.time(); 
show(toc-tic); 
# 0.09 seconds 

簡單的檢查，顯示結果是相同的：

all(filter == filter2) 
# TRUE