Что вызывает разницу между cal c и cellStats в растровых вычислениях в R? - PullRequest
0 голосов
/ 19 апреля 2020

Я работаю с набором данных, который состоит из 20 слоев, уложенных в RasterBrick (происходящий из массива). Я посмотрел на сумму слоев, рассчитанных как с «Cal c» и «CellStats». Я использовал cal c для вычисления суммы общих значений и cellStats для просмотра среднего значения на слой (полезно для временного ряда). Однако, когда я суммирую среднее значение каждого слоя, оно вдвое меньше значения другой вычисленной суммы. Что вызывает эту разницу? Что я пропускаю?

Код выглядит следующим образом:

testarray <- runif(54214776,0,1) 
# Although testarray should contain a raster of 127x147 with 2904 time layers. 
# Not really sure how to create that yet. 

for (i in 1830:1849){
  slice<-array2[,,i]
  r <- raster(nrow=(127*5), ncol=(147*5), resolution =5, ext=ext1, vals=slice)
  x <- stack(x , r)
}

brickhp2 <- brick(x)

r_sumhp2 <- calc(brickhp2, sum, na.rm=TRUE)
r_sumhp2[r_sumhp2<= 0] <- NA


SWEavgpertimestepM <- cellStats(brickhp2, stat='mean', na.rm=TRUE)

Цель состоит в том, чтобы сравнить сумму слоев, рассчитанную с помощью 'cal c (x, sum)', с сумма среднего значения, рассчитанная с помощью cellStats (x, mean).

Rasterbrick выглядит следующим образом (600kb, GTiff): http://www.filedropper.com/brickhp2 * Если есть лучший способ поделиться этим Пожалуйста, дайте мне знать.

1 Ответ

1 голос
/ 05 мая 2020

Путаница возникает, когда вы используете calc, который работает по пикселям для кирпича (т.е. выполняет вычисление по 20 значениям в каждом пикселе и возвращает один растровый слой) и cellStats, который выполняет вычисления для каждого растровый слой индивидуально и возвращает отдельные значения для каждого слоя. Вы можете увидеть, что результаты сопоставимы, если вы используете этот код:

library(raster)
##set seed so you get the same runif vals
set.seed(999)

##create example rasters
ls=list()
for (i in 1:20){
  r <- raster(nrow=(127*5), ncol=(147*5), vals=runif(127*5*147*5))
  ls[[i]] <- r
}
##create raster brick
brickhp2 <- brick(ls)

##calc sum (pixel-wise)
r_sumhp2 <- calc(brickhp2, sum, na.rm=TRUE)
r_sumhp2 ##returns raster layer

##calc mean (layer-wise)
r_meanhp2 <- cellStats(brickhp2, stat='mean', na.rm=TRUE)
r_meanhp2 ##returns vector of length nlayers(brickhp2)

##to get equivalent values you need to divide r_sumhp2 by the number of layers 
##and then calculate the mean
cellStats(r_sumhp2/nlayers(brickhp2),stat="mean")

[1] 0,4999381

##and for r_meanhp2 you need to calculate the mean of the means
mean(r_meanhp2)

[1] 0,4999381

Вам потребуется определите для себя, хотите ли вы использовать пиксельный или послойный результат для своего приложения.

...