опция agstudy - это здорово. Другой вариант заключается в том, чтобы заранее создать карту итераций, а затем применить применительно к этому:
iterate <- expand.grid(1:3,1:3)
X <- apply(iterate,1,function(x) myfunction(ddf[,x[1]],ddf[,x[2]]))
Возможно, стоит спуститься по этому маршруту, только если у вас имеется большое количество столбцов и многоядерная система. Затем вы будете параллелизировать его, используя одну из операций параллельного применения в пакете parallel
.