# Speed up a sampling function

I have the following dataset

kkk<-data.frame(days=1:100,positive=rbinom(100,1,0.05))


Over the monitoring period of 100 days, if an event occurs then for that day kkk$positive==1 else kkk$positive==0

The following function, samples k=s.freq=1:12 number of times j=s.length=1:30 consecutive elements from the sample space and checks if at least in one of the samples an event was recorded. The sample space consists of all possible j consecutive combinations of kkk$days. The elements of the sample space are overlapping. Before each of the k=s.freq=1:12 samples, the sample space is reevaluated and elements that have common days with the sampled element are removed before the next sampling occurs (for the total of k samplings). my.sampling<-function(days, status, s.length, s.freq, iter){ start<-Sys.time() #build the data frame ddd<-data.frame(days=days, positive=status) #length of each sample j<-s.length #Buld the sample space monitor.ss<-matrix(nrow=dim(ddd)[1]-j+1, ncol=j+1) n.row<-dim(ddd)[1]-j+1 #Create the elements of sample space: each row of the monitor.ss is an element of the sample space for (i in 1:j){ monitor.ss[,i]<-i:(n.row-1+i) } rm(i) #adds if the each of the possible samples of the sample space at least one event was observed or not monitor.ss[,j+1]<-apply(monitor.ss,1,function(x) { r.low<-range(x,na.rm=T)[1] r.high<-range(x,na.rm=T)[2] return(as.numeric(sum(ddd$positive[r.low:r.high])>0))})

#Build the initial sample space as data frame
sampling.start<-data.frame(monitor.ss)

suc.list<-rep(NA, iter) #Initiate a vector to keep track of successes

#Now the sampling takes place

for (i in 1:iter) { #for each iteration

sampling.start<-data.frame(monitor.ss)
k<-s.freq #number of samples drawn from the sample space
k.t=0 #controller to break the while loop
suc<-0

while (k.t<=k | dim(sampling.start)[1]>0) { #breaks the while loop if k reaches the specified limit (k=s.freq) or if the reduction of the sample space (see below) does not leave any elements in the sample space

s<-sampling.start[sample(nrow(sampling.start),1),] #samples the first trial
suc<-suc+s[1,j+1] #adds up if there is a success or not ( a tleast one event in the sample was observed)
sampling.start <- sampling.start[apply(sampling.start, 1, function(x) !any(x %in% as.numeric(s[,1:j]))),] #sample space reduction: Elements of the sample space that intersect with the sample that was just selected are removed from the sample space before the next sampling occurs
k.t<-k.t+1 #controller for the wile loop
}

suc.list[i]<-(suc/k.t)

}
print(Sys.time()-start)
message("")
#hist(suc.list, breaks=sqrt(iter))
return(suc.list)
}


For example:

my.sampling(kkk$days, kkk$positive, s.length=3, s.freq=5, iter=300)


My problem is that the function requires about 100msec for each iteration and I need to run this function about 100.000 times with at least 1000 iterations each time. I have used preallocation of the vectors and apply to speed things up, however its still too slow.

Parallel processing is applicable to my problem as I can split the data on which the function will have to run on multiple cores, however if my calculations are right I need about 120 days of processing pro CPU core.

Is this possible to reduce the run time of the function? I have tried cpmfun but without much improvement. Otherwise I have to resort to some type of AWS EC2 solution.

Thanks a lot

• You can, for one, replace that entire allocation of monitor.ss where you are "building sample space" with this single line: monitor.ss <- embed(1:dim(ddd)[1], j)[ , j:1]. I haven't been able to follow your problem description, however.
– DWin
Commented Jun 21, 2012 at 17:24
• Nice @DWin, However that loop doen't take any time. His main problen is in the last loop. I couldn't follow either.
– lselzer
Commented Jun 21, 2012 at 17:38
• Thanks @DWin. However as Iselzer writes this parts does not consume significant CPU time. I also think the last part is the most CPU intensive. The most important part of the sampling is the reevaluation and reduction of the sample space (the while loop) before the next consecutive sampling. Any ideas if the function or this part can get faster?
– ECII
Commented Jun 21, 2012 at 22:26

This question is probably not appropriate for stackoverflow, since it requires mainly code review.

First of all, there is a logics problem with your while loop. Changing | to & reduces the number of loops considerably.

Furthermore you should not use dataframes, if you can use a matrix.

my.sampling1<-function(days, status, s.length, s.freq, iter){
start<-Sys.time()

#build the data frame
ddd<-data.frame(days=days, positive=status)

#length of each sample
j<-s.length

#Buld the sample space
monitor.ss<-matrix(nrow=dim(ddd)[1]-j+1, ncol=j+1)
n.row<-dim(ddd)[1]-j+1

#Create the elements of sample space: each row of the monitor.ss is an element of the sample space
for (i in 1:j){
monitor.ss[,i]<-i:(n.row-1+i)
}
rm(i)

#adds if the each of the possible samples of the sample space at least one event was observed or not
monitor.ss[,j+1]<-apply(monitor.ss,1,function(x) {
r.low<-range(x,na.rm=T)[1]
r.high<-range(x,na.rm=T)[2]
return(as.numeric(sum(ddd$positive[r.low:r.high])>0))}) suc.list<-rep(NA, iter) #Initiate a vector to keep track of successes #Now the sampling takes place k<-s.freq #number of samples drawn from the sample space for (i in 1:iter) { #for each iteration #do not convert to data.frame sampling.start<-monitor.ss k.t=0 #controller to break the while loop suc<-0 #use & to stop the loop if one condition is not fullfilled, second condition changed because a matrix is used instead of a data.frame while (k.t<=k & !is.null(nrow(sampling.start))) { #breaks the while loop if k reaches the specified limit (k=s.freq) or if the reduction of the sample space (see below) does not leave any elements in the sample space #set.seed(k.t*i) #use set.seed to test if results are identical s<-sampling.start[sample(nrow(sampling.start),1),] #samples the first trial suc<-suc+s[j+1] #adds up if there is a success or not ( a tleast one event in the sample was observed) sampling.start <- sampling.start[apply(sampling.start, 1, function(x) !any(x %in% as.numeric(s[1:j]))),] #sample space reduction: Elements of the sample space that intersect with the sample that was just selected are removed from the sample space before the next sampling occurs k.t<-k.t+1 #controller for the wile loop } suc.list[i]<-(suc/k.t) } print(Sys.time()-start) message("") #hist(suc.list, breaks=sqrt(iter)) return(suc.list) } kkk<-data.frame(days=1:100,positive=rbinom(100,1,0.05)) res1<-my.sampling1(kkk$days, kkk\$positive, s.length=3, s.freq=5, iter=20)


This is already faster. To make it even more efficient, I would try to avoid apply inside the loop (could an algebraic operation be used instead?) and try to vectorize. Furthermore you can tidy up. There are some unnecessary assignments inside the function.