1

シェープファイルの各ポリゴンに対して計算された平均 +/- 2.5 std を超える値を削除して、ラスター ファイルを標準化しようとしています。ゾーン統計を計算してシェープファイルにフィードバックすることはできましたが、シェープファイルに含まれるデータでラスター ファイルをマスクすることができません。コードの最後の 2 行を実行すると、エラー メッセージが表示されます。平均 + 2.5 std を超える値を抽出すると、出力としてラスターではなく値 (NA) のみが得られ、その結果、マスク関数は次のエラーを返します。

(関数 (classes、fdef、mtable) のエラー: シグネチャ '"RasterBrick"、"logical"' の関数 'mask' の継承されたメソッドが見つかりません</p> 

#Read raster file
library(raster)
prec <- getData('worldclim', var='prec', res=10)

#Read shapefile
france <- getData('GADM', country='FRA', level=1)

#Convert shapefile to raster
r_fra <- rasterize(france, prec)

#Calculate mean for each polygon
mean_fra_prec <- zonal(prec, r_fra, FUN='mean')

#Calculate standard deviation for each polygon
std_fra_prec <- zonal(prec, r_fra, FUN = 'sd')

#Calculate mean +/- 2.5 std for each polygon
ul_prec  <- mean_fra_prec + (2.5 * std_fra_prec)
ll_prec <- mean_fra_prec - (2.5*std_fra_prec)

#Merge mean +/- 2.5 std for each polygon
x <- data.frame(zone=1:22)
ul_prec2 <- merge(x, ul_prec, by='zone', all.x=TRUE)
ll_prec2 <- merge(x, ul_prec, by='zone', all.x=TRUE)

#Join new data to shapefile
france@data <- cbind(france@data, ul_prec2[,-1])
france@data <- cbind(france@data, ll_prec2[,-1])

#Remove values above mean + 2.5 std for each polygon
ext_prec <- prec[!(prec > france@data$ul_prec2)] <- NA #Values above mean + 2.5 std
ma_prec <- mask(prec, ext_prec) #Mask values above mean + 2.5 std
4

1 に答える 1

2

ここに作業スクリプトがあると思います。

library(raster)
# Get raster data
prec <- getData('worldclim', var='prec', res=10)
# Get polygons
france <- getData('GADM', country='FRA', level=1)

# for this example, simplify to one layer, smaller extent
prec <- crop(prec[[1]], france)

# rasterize polygons
r_fra <- rasterize(france, prec)

# Calculate mean for each polygon
mean_fra_prec <- zonal(prec, r_fra, fun='mean')
# order and remove ID
mean_fra_prec <- mean_fra_prec[order(mean_fra_prec[,1]), -1]
#Calculate standard deviation for each polygon
std_fra_prec <- zonal(prec, r_fra, fun = 'sd')
std_fra_prec <- std_fra_prec[order(std_fra_prec[,1]), -1]

# alternative route, may fail on large data set
# v <- extract(prec, france)
# mean_fra_prec <- sapply(v, mean)
# std_fra_prec <- sapply(v, sd)

#Calculate mean +/- 2 std for each polygon
ul  <- mean_fra_prec + 2*std_fra_prec
ll <- mean_fra_prec - 2*std_fra_prec
d <- data.frame(1:length(ll), ll, ul)

v <- subs(r_fra, d, which=2:3)

keep <- prec > v[[1]] & prec < v[[2]]

x <- mask(prec, keep, maskvalue=FALSE)
plot(x)
于 2013-06-13T07:54:26.200 に答える