Это должно сработать в Base-R, у кого-то еще может быть «более чистое» решение из пакета.
DFSpecies <- cbind(DFSpecies[1:4],NA,NA,NA,NA,NA,NA,NA,NA,DFSpecies[5:80])
for(i in 1:19){
DFSpecies <- cbind(DFSpecies[1:((i*12)+4)],NA,NA,NA,NA,NA,NA,NA,NA,DFSpecies[(i*12+1):NCOL(DFSpecies)])
}
DFSpecies <- cbind(DFSpecies,NA,NA,NA,NA,NA,NA,NA,NA)