减少大量栅格上马赛克的内存使用
Reduce memory usage for mosaic on large list of rasters
我正在使用 raster
包中的 mosaic
函数,使用@RobertH here 建议的方法组合一长串(11,000 个文件)光栅。
rlist <- sapply(list_names)
rlist$fun <- mean
rlist$na.rm <- TRUE
x <- do.call(mosaic, rlist)
如您所想,这最终会超出我的可用内存(在几个不同的机器和计算集群上)。我的问题是:有没有办法减少 mosaic
或 do.call
的内存使用量?我试过在 rasterOptions()
中更改 maxmemory
,但这似乎没有帮助。以较小的批次处理栅格似乎有问题,因为栅格可能在空间上是分离的(即,顺序栅格文件可能彼此相距很远)。提前感谢您提供的任何帮助。
不是一次将所有栅格加载到内存中(在 mosaic()
调用中),您可以一次处理一个吗?这样,每次将更多光栅存入内存时,您的马赛克都会更新,但随后您可以摆脱新光栅,只保留不断更新的马赛克光栅。
假设您的 rlist
对象是栅格列表,我在想类似的东西:
伪代码
- 将一个
updating_raster
对象初始化为列表中的第一个栅格
- 从第二个光栅开始,依次循环遍历列表中的每个光栅
- 将第 i 个栅格读入内存
next_raster
- 更新
updating_raster
对象,方法是用其自身的马赛克和使用加权平均值的下一个光栅覆盖它
R
代码
正在使用 mosaic()
帮助文件示例中的代码进行测试...
首先生成一些栅格并使用标准的镶嵌方法。
library(raster)
r <- raster(ncol=100, nrow=100)
r1 <- crop(r, extent(-10, 11, -10, 11))
r2 <- crop(r, extent(0, 20, 0, 20))
r3 <- crop(r, extent(9, 30, 9, 30))
r1[] <- 1:ncell(r1)
r2[] <- 1:ncell(r2)
r3[] <- 1:ncell(r3)
m1 <- mosaic(r1, r2, r3, fun=mean)
将栅格放入列表中,使它们的格式与我认为的相似。
rlist <- list(r1, r2, r3)
由于 NA
处理 weighted.mean()
函数,我选择通过将求和和除法分解为不同的步骤来创建相同的效果...
首先初始化求和栅格:
updating_sum_raster <- rlist[[1]]
然后初始化"counter"光栅。这将代表在每个像素处进入镶嵌的栅格数量。在不是 NA
的所有单元格中,它都以 1 开头。它应该正确处理 NA
s,这样它只会在给定像素上增加一个非 NA
值被添加到更新总和。
updating_counter_raster <- updating_sum_raster
updating_counter_raster[!is.na(updating_counter_raster)] <- 1
这是一个不需要所有栅格一次都在内存中的循环。添加到马赛克的栅格的计数器栅格的值仅在不是 NA
的像元中为 1。通过将当前计数器栅格和更新计数器栅格相加来更新计数器。通过对当前栅格值和更新栅格值求和来更新总和。
for (i in 2:length(rlist)) {
next_sum_raster <- rlist[[i]]
next_counter_raster <- next_sum_raster
next_counter_raster[!is.na(next_counter_raster)] <- 1
updating_sum_raster <- mosaic(x = updating_sum_raster, y = next_sum_raster, fun = sum)
updating_counter_raster <- mosaic(updating_counter_raster, next_counter_raster, fun = sum)
}
m2 <- updating_sum_raster / updating_counter_raster
这里的值似乎与mosaic()
函数的使用相匹配
identical(values(m1), values(m2))
> TRUE
但栅格本身并不相同:
identical(m1, m2)
> FALSE
不太清楚为什么,但也许这会让你更接近?
也许 compareRaster()
是更好的检查方法:
compareRaster(m1, m2)
> TRUE
万岁!
剧情来了!
plot(m1)
text(m1, digits = 2)
plot(m2)
text(m2, digits = 2)
在杂草中多挖一点...
来自 mosaic.R 文件:
看起来 mosaic()
函数初始化了一个名为 v
的矩阵,以填充列表中所有栅格中所有像元的值。矩阵中的行数 v
是输出栅格中的像元数(基于完整镶嵌范围和分辨率),列数是您的情况下要镶嵌的栅格数(11,000) .也许您 运行 进入了 R 中矩阵创建的极限?
对于 1000 x 1000 光栅(1e6 像素),NA
的 v
矩阵占用 41 GB。您希望最终的镶嵌栅格有多大?
r <- raster(ncol=1e3, nrow=1e3)
x <- 11000
v <- matrix(NA, nrow=ncell(r), ncol=x)
format(object.size(v), units = "GB")
[1] "41 Gb"
我正在使用 raster
包中的 mosaic
函数,使用@RobertH here 建议的方法组合一长串(11,000 个文件)光栅。
rlist <- sapply(list_names)
rlist$fun <- mean
rlist$na.rm <- TRUE
x <- do.call(mosaic, rlist)
如您所想,这最终会超出我的可用内存(在几个不同的机器和计算集群上)。我的问题是:有没有办法减少 mosaic
或 do.call
的内存使用量?我试过在 rasterOptions()
中更改 maxmemory
,但这似乎没有帮助。以较小的批次处理栅格似乎有问题,因为栅格可能在空间上是分离的(即,顺序栅格文件可能彼此相距很远)。提前感谢您提供的任何帮助。
不是一次将所有栅格加载到内存中(在 mosaic()
调用中),您可以一次处理一个吗?这样,每次将更多光栅存入内存时,您的马赛克都会更新,但随后您可以摆脱新光栅,只保留不断更新的马赛克光栅。
假设您的 rlist
对象是栅格列表,我在想类似的东西:
伪代码
- 将一个
updating_raster
对象初始化为列表中的第一个栅格 - 从第二个光栅开始,依次循环遍历列表中的每个光栅
- 将第 i 个栅格读入内存
next_raster
- 更新
updating_raster
对象,方法是用其自身的马赛克和使用加权平均值的下一个光栅覆盖它
R
代码
正在使用 mosaic()
帮助文件示例中的代码进行测试...
首先生成一些栅格并使用标准的镶嵌方法。
library(raster)
r <- raster(ncol=100, nrow=100)
r1 <- crop(r, extent(-10, 11, -10, 11))
r2 <- crop(r, extent(0, 20, 0, 20))
r3 <- crop(r, extent(9, 30, 9, 30))
r1[] <- 1:ncell(r1)
r2[] <- 1:ncell(r2)
r3[] <- 1:ncell(r3)
m1 <- mosaic(r1, r2, r3, fun=mean)
将栅格放入列表中,使它们的格式与我认为的相似。
rlist <- list(r1, r2, r3)
由于 NA
处理 weighted.mean()
函数,我选择通过将求和和除法分解为不同的步骤来创建相同的效果...
首先初始化求和栅格:
updating_sum_raster <- rlist[[1]]
然后初始化"counter"光栅。这将代表在每个像素处进入镶嵌的栅格数量。在不是 NA
的所有单元格中,它都以 1 开头。它应该正确处理 NA
s,这样它只会在给定像素上增加一个非 NA
值被添加到更新总和。
updating_counter_raster <- updating_sum_raster
updating_counter_raster[!is.na(updating_counter_raster)] <- 1
这是一个不需要所有栅格一次都在内存中的循环。添加到马赛克的栅格的计数器栅格的值仅在不是 NA
的像元中为 1。通过将当前计数器栅格和更新计数器栅格相加来更新计数器。通过对当前栅格值和更新栅格值求和来更新总和。
for (i in 2:length(rlist)) {
next_sum_raster <- rlist[[i]]
next_counter_raster <- next_sum_raster
next_counter_raster[!is.na(next_counter_raster)] <- 1
updating_sum_raster <- mosaic(x = updating_sum_raster, y = next_sum_raster, fun = sum)
updating_counter_raster <- mosaic(updating_counter_raster, next_counter_raster, fun = sum)
}
m2 <- updating_sum_raster / updating_counter_raster
这里的值似乎与mosaic()
函数的使用相匹配
identical(values(m1), values(m2))
> TRUE
但栅格本身并不相同:
identical(m1, m2)
> FALSE
不太清楚为什么,但也许这会让你更接近?
也许 compareRaster()
是更好的检查方法:
compareRaster(m1, m2)
> TRUE
万岁!
剧情来了!
plot(m1)
text(m1, digits = 2)
plot(m2)
text(m2, digits = 2)
在杂草中多挖一点...
来自 mosaic.R 文件:
看起来 mosaic()
函数初始化了一个名为 v
的矩阵,以填充列表中所有栅格中所有像元的值。矩阵中的行数 v
是输出栅格中的像元数(基于完整镶嵌范围和分辨率),列数是您的情况下要镶嵌的栅格数(11,000) .也许您 运行 进入了 R 中矩阵创建的极限?
对于 1000 x 1000 光栅(1e6 像素),NA
的 v
矩阵占用 41 GB。您希望最终的镶嵌栅格有多大?
r <- raster(ncol=1e3, nrow=1e3)
x <- 11000
v <- matrix(NA, nrow=ncell(r), ncol=x)
format(object.size(v), units = "GB")
[1] "41 Gb"