Более быстрые способы создания матриц в R

Question

У меня есть следующий пример того, что я делал, что формально просто, но я хотел проверить, каковы потенциальные альтернативы моему коду - чтобы ускорить, если это возможно. Вот пример:

Time1=Sys.time()
v=rep(c("A","B","C","D","E","F","G","H","I","J","K","L","M","N","O","P"),
     each=1000)
m=matrix(0,ncol=length(v),nrow=length(v))
for (j in 1:length(v)) {
  for(i in 1:length(v)) {
    if (v[j]==v[i]) {
      m[j,i]=1
    } else {
      next
    }
  }
}
Time2=Sys.time()
Time2-Time1
# Time difference of 1.405404 mins

Я создаю простую реляционную матрицу, где вектор v1 можно интерпретировать как помещенный в виде строк и столбцов, а в матричные карты, где результаты равны. Если они равны, мы получаем m[j,i]=1; если не равно, m[j,i]=0. Как я уже сказал, я хотел бы сделать этот код go быстрее. Я пытался придумать способ кодировать его как apply функцию, но пока не понял этого. Тем не менее, я хотел бы знать, есть ли другие варианты помимо того, что я сказал.

РЕДАКТИРОВАТЬ : я внес в текст некоторые исправления и попытался уточнить вопрос.

Answer 1

Конечно. Если предположить, что ваши выборочные данные не являются репрезентативными для реальных данных, то это работает примерно в 6 раз быстрее:

m2 <- +(outer(v, v, `==`))
all.equal(m, m2)
# [1] TRUE

Если, однако, ваши реальные данные имеют значительное количество дубликатов, то @ Sathi sh ' s метод удаления дубликатов перед сравнением и распространения их через матрицу скорее всего намного быстрее.

Answer 2

Я брошу свою шляпу в кольцо с CJ из data.table.

libary(data.table)
m2 <- matrix(+(CJ(v1 = v,v2 = v,sorted=FALSE)[,ans := v1==v2][,ans]),length(v))
all.equal(m,m2)
#[1] TRUE

Answer 3

Я думаю, outer подход @ r2evans - самый простой способ построения матрицы. Ниже приведена еще одна базовая опция R с использованием expand.grid

m2 <- matrix(+do.call("==",expand.grid(v,v)),length(v))

Answer 4

Если у вас много нулей, вы можете использовать sparseMatrix. Вы заполняете позиции, которые совпадают, и оставляете остальные нули. Это потребляет меньше памяти, но вы можете использовать его только с определенными функциями. Он используется, например, в gl mnet lasso.

Я думаю, что решение @ r2evans является наиболее лаконичным и будет работать в большинстве случаев.

Ниже у меня есть несколько кодов из ответов, большинство из них действительно быстрее, чем ОП

library(microbenchmark)
library(Matrix)
library(data.table)
setDTthreads(threads =1)

f_sw = function(v){
N = length(v)
i = lapply(v,function(i)which(v==i))
j = rep(1:N,times=sapply(i,length))
as.matrix(sparseMatrix(i=unlist(i),j=j,dims=list(N,N)))
}

f_r2evans = function(v){
  m2 <- +(outer(v, v, `==`))
  return(m2)
}

f_IanCampbell = function(v){ 
  matrix(+(CJ(v1 = v,v2 = v,sorted=FALSE)[,ans := v1==v2][,ans]),length(v))
}

microbenchmark(f_IanCampbell(v),f_sw(v),f_r2evans(v),times=5)

Unit: seconds
             expr       min        lq      mean    median        uq      max
 f_IanCampbell(v) 10.820616 11.325422 12.544146 12.983926 13.126655 14.46411
          f_sw(v)  7.014364  7.228585  8.206858  7.745741  8.877425 10.16818
     f_r2evans(v)  9.117405  9.519443  9.996789  9.896823 10.288586 11.16169
 neval cld
     5   b
     5  a 
     5  a

Answer 5

может быть таким:

Вместо проверки на равенство всех возможных значений, мы можем сделать это для вектора и повторить матрицу из него до 1000 раз по строкам и столбцам. Это даст тот же результат. Порядок столбцов и строк не поддерживается этим кодом. Но, используя имена строк и столбцов, мы можем убедиться, что ответ правильный или нет.

Я использовал t(), потому что привязка столбцов происходит быстрее, чем привязка строк.

system.time({
  v1 <- c("A","B","C","D","E","F","G","H","I","J","K","L","M","N","O","P")
  m1 <- sapply(v1, function(x) as.integer(v1 == x))
  rownames(m1) <- colnames(m1)
  m1 <- do.call('cbind', mget(rep('m1', 1000)))
  m1 <- t(m1)
  m1 <- do.call('cbind', mget(rep('m1', 1000)))
  m1 <- t(m1)
}) 
# user  system elapsed 
# 9.32    0.50    9.84 


dim(m1)
# [1] 16000 16000

---

Другой метод:

Этот (жестко запрограммированный) не будет выполнять никакого сравнения, но мы создаем значения на основе того, что может произойти, сравнивая вектор с его собственными значениями.

Его можно улучшить с помощью конструкции eval(parse(text=paste())).

system.time({
  m4 <- matrix(data = 
           c(
             c(rep(1, 1000), rep(0, 15000)),
             c(rep(0, 1000), rep(1, 1000), rep(0, 14000)),
             c(rep(0, 2000), rep(1, 1000), rep(0, 13000)),
             c(rep(0, 3000), rep(1, 1000), rep(0, 12000)),
             c(rep(0, 4000), rep(1, 1000), rep(0, 11000)),
             c(rep(0, 5000), rep(1, 1000), rep(0, 10000)),
             c(rep(0, 6000), rep(1, 1000), rep(0, 9000)),
             c(rep(0, 7000), rep(1, 1000), rep(0, 8000)),
             c(rep(0, 8000), rep(1, 1000), rep(0, 7000)),
             c(rep(0, 9000), rep(1, 1000), rep(0, 6000)),
             c(rep(0, 10000), rep(1, 1000), rep(0, 5000)),
             c(rep(0, 11000), rep(1, 1000), rep(0, 4000)),
             c(rep(0, 12000), rep(1, 1000), rep(0, 3000)),
             c(rep(0, 13000), rep(1, 1000), rep(0, 2000)),
             c(rep(0, 14000), rep(1, 1000), rep(0, 1000)),
             c(rep(0, 15000), rep(1, 1000))), nrow = 16000, ncol = 16000)
})
# user  system elapsed 
# 0.72    0.93    1.82 

Примечание: Как сказал @ r2evans, это не будет работать, если выборочные данные ОП не являются репрезентативными для реальных данных