0

我想知道是否有更有效的方法来加速下面的代码。它使用嵌套循环中需要子采样的过程（先前的答案https://stackoverflow.com/a/13629611/1176697有助于提高效率）。当 B=500 时，R 倾向于挂起，尽管计算机操作系统并未受到过度影响。

目标是使用 B=1000 并使用更大的 m 值（m=75,m=100,m=150）运行以下代码

我在下面的代码中详细说明了该过程，并包含了指向可重现数据集的链接。

#Estimation of order m `leave one out' hyperbolic efficiency scores

#The procedure sequentially works through each observation in `IOs' and 
#calculates a DEA order M efficiency score by leaving out the observation
# under investigation in the DEA reference set

# Step 1: Load the packages, create Inputs (x) and Outputs (y), choose 
# m(the order of the partial frontier or reference set to be used in DEA)
# and B the number of monte carlo simulations of each m order DEA estimate

# Step 2: For each observations a in x, x1 and y1
#are create which 'leaves out' this observation. 

# Step 3: From these matrices subsamples (xref, yref) of size [m,] are 
# taken and used in DEA estimation.

# Step 4: The DEA estimation uses the m subsample from step 3
# as a reference set and evaluates the efficiency of the observation that
# has been 'left out'
#(thus the first two arguments in DEA are matrices of order [1,3] )

# Step 5: Steps 3 and 4 are repeated B times to obtain B simulations of the 
# order m efficiency score and a mean and standard deviation are
# calculated and placed in effm.

# IOs data can be found here: https://dl.dropbox.com/u/1972975/IOs.txt
# From IOs an Input matrix (x[1376,3]) and an Output matrix (y[1376,3])
# are created.

library(Benchmarking)
x <- IOs[,1:3]
y<-IOs[,4:6]
A<-nrow(x)
effm <- matrix(nrow = A, ncol = 2)
m <- 50
B <- 500
pb <- txtProgressBar(min = 0,
                     max = A, style=3)

for(a in 1:A) {
  x1 <- x[-a,]
  y1 <- y[-a,]
  theta <- numeric(B)
  xynrow<-nrow(x1)
  mB<-m*B 
  xrefm <- x1[sample(1:xynrow, mB, replace=TRUE),] # get all of your samples at once(https://stackoverflow.com/a/13629611/1176697)
  yrefm <- y1[sample(1:xynrow, mB, replace=TRUE),]
  deaX <- as.matrix(x[a,], ncol=3)
  deaY <-as.matrix(y[a,], ncol=3)

  for(i in 1:B){
    theta[i] <- dea(deaX, deaY, RTS = 'vrs', ORIENTATION = 'graph',
                    xrefm[(1:m) + (i-1) * m,], yrefm[(1:m) + (i-1) * m,], FAST=TRUE)
  }

  effm[a,1] <- mean(theta)
  effm[a,2] <- sd(theta) / sqrt(B)
  setTxtProgressBar(pb, a) 
}
close(pb)
4

0 回答 0