处理大数据集计算R中的SPI
我正在处理大型全球降水数据集（具有非常精细的空间分辨率），以便使用R中的SPEI包计算标准化降水指数。
总的来说，我的问题是指使用非常大的数据优化数据处理。 我在其他帖子中发现了一些讨论（这里，在这里，这里是fo实例），但没有什么与我的情况类似。
我的输入是一个矩阵，其中超过20年的降水时间序列的月观测值（> 20 * 12行）为> 1,000,000点（列）。 SPI的计算为每个时间序列执行一系列步骤，并将该指数计算为中值的标准偏差。 输出是一个列表，结果矩阵（$拟合）具有相同的输入矩阵大小。
这里是一个代码示例：
require(SPEI)
#generating a random values matrix
data<-replicate(200, rnorm(240))
# erasing negative values
data[data<=0]=0
spi6 <- spi(data, 6, kernel = list(type = 'rectangular', shift = 0), distribution = 'PearsonIII', fit = 'ub-pwm', na.rm = FALSE, ref.start=NULL, ref.end=NULL, x=FALSE, params=NULL)
#testing the results
plot(spi6$fitted[,67])
#taking my results out
results <- t(spi6$fitted)
这个脚本完美地工作，但是如果我增加点数（在这种情况下为列），处理时间将呈指数增长。 直到遇到内存不足问题：
Warning messages:
1: In std[ff, s] <- qnorm(cdfpe3(acu.pred[ff], p3par)) :
Reached total allocation of 16253Mb: see help(memory.size)
2: In std[ff, s] <- qnorm(cdfpe3(acu.pred[ff], p3par)) :
Reached total allocation of 16253Mb: see help(memory.size)
3: In NextMethod("[<-") :
Reached total allocation of 16253Mb: see help(memory.size)
4: In NextMethod("[<-") :
Reached total allocation of 16253Mb: see help(memory.size)
5: In std[ff, s] <- qnorm(pze + (1 - pze) * pnorm(std[ff, s])) :
Reached total allocation of 16253Mb: see help(memory.size)
6: In std[ff, s] <- qnorm(pze + (1 - pze) * pnorm(std[ff, s])) :
Reached total allocation of 16253Mb: see help(memory.size)
7: In NextMethod("[<-") :
Reached total allocation of 16253Mb: see help(memory.size)
8: In NextMethod("[<-") :
Reached total allocation of 16253Mb: see help(memory.size)
如何拆分我的输入矩阵（或者将输入矩阵拆分过程），以便并行处理列向量组（每个列向量都是特定点的全部时间序列），而不会丢失信息（或者绕过我的数据）？ 谢谢。
我发现了一个非常线性的，而不是指数级的处理时间：
system.time(spi(data[,1], 6))
system.time(spi(data[,1:10], 6))
system.time(spi(data[,1:100], 6))
如果您看到指数级增长，则可能是由于RAM分配过多的问题导致的。
一个简单的解决方案是在矩阵上分割计算：
spi6 <- data*NA
system.time(
for (i in 1:100) spi6[,i] <- spi(data[,i], 6)$fitted
)
或者具有相似的效率：
system.time(
spi6 <- apply(data[,1:100], 2, function(x) spi(x, 6)$fitted)
)
正如您所看到的，计算时间与原始选项非常相似，您可以将整个矩阵作为spi()函数的输入。 但是 ，如果您遇到内存问题，这可能会解决它们。
另一方面，如果您有权访问多核计算机（因为现在很可能是这种情况），您可以通过并行化计算来看到计算时间的改进：
library(snowfall)
sfInit(parallel=TRUE, cpus=2)
sfLibrary(SPEI)
system.time(
spi6 <- sfApply(data[,1:100], 2, function(x) spi(x, 6)$fitted)
)
sfStop()
您可能需要将更高的值设置为ncpu以获得更高的速度增益，具体取决于您的计算机支持的线程数。 但是 ， sfApply并不能解决大型数据集的内存问题。 这是因为函数会将数据集分配给分配的cpus数量。 由于系统的总内存不会改变，因此您的内存不足也会相同。
解决方案是合并两种方法：拆分数据集，然后并行化。 像这样的东西：
data <- replicate(10000, rnorm(240))
sfInit(parallel=TRUE, cpus=10)
sfLibrary(SPEI)
spi6 <- data*NA
for (i in 0:9) {
chunk <- (i*1000)+(1:1000)
spi6[,chunk] <- sfApply(data[,chunk], 2, function(x) spi(x, 6)$fitted)
}
sfStop()
现在，您只需要查找块的最大大小，即传递给sfApply的数据sfApply ，以便不会溢出可用RAM。 有一些尝试和错误很容易。