首页
学习
活动
专区
圈层
工具
发布
首页
学习
活动
专区
圈层
工具
MCP广场
社区首页 >问答首页 >使用多个函数从栅格中提取值

使用多个函数从栅格中提取值
EN

Stack Overflow用户
提问于 2021-02-15 17:02:32
回答 2查看 31关注 0票数 0

我在这篇文章中使用了解决方案中的代码:

How to pass multiple function in raster::extract()

代码语言:javascript
运行
复制
functions_summ <- function(x, na.rm) c(mean = mean(x, na.rm=na.rm), min = min(x, na.rm=na.rm), max = max(x, na.rm=na.rm), sd = sd(x, na.rm = na.rm))

extraction <- extract( Stack_for_Classification,TrVal, fun=functions_summ, na.rm=TRUE, df = TRUE, sp = T) ## the sp = T. means that I am adding polygon values to the csv

nom <- sapply(TrVal@polygons, slot, "ID")

extraction <- data.frame(ID = nom, Value = extraction)

names(extraction) <- gsub(x = names(extraction), pattern = "\\Value.", replacement = "")  
typeof(extraction)

View(extraction)

Stack_for_Classification是堆叠栅格,其中包含我想要使用多边形提取的值。TrVal是包含多边形的shapefile。它有一个带有三个标识值的字段"Type“,我用它来根据提取的值来查看三种不同类型之间的差异。问题是,在前面的代码中,我得到了均值、最大值、最小值和后缀为"0,1,2“的sd。总共12。是否可以使用字段“类型”来命名行,例如"mean_Type1","max_Type2"..诸若此类?使用View(提取)后上传图片

View(extraction) result

EN

回答 2

Stack Overflow用户

发布于 2021-04-23 21:36:23

你可以这样做:

代码语言:javascript
运行
复制
library(raster)

# simulate data -------------------------------------------------------------
r <- s <- t <- raster()
r[] <- 1:ncell(r)
s[] <- sample(1:10,ncell(s),replace = T)
t[] <- runif(ncell(t))

stacked <- stack(r,s,t)

#polys (from the SpatialPolygonsDataFrame example)
Sr1 = Polygon(cbind(c(2,4,4,1,2),c(2,3,5,4,2)))
Sr2 = Polygon(cbind(c(5,4,2,5),c(2,3,2,2)))
Sr3 = Polygon(cbind(c(4,4,5,10,4),c(5,3,2,5,5)))
Sr4 = Polygon(cbind(c(5,6,6,5,5),c(4,4,3,3,4)), hole = TRUE)

Srs1 = Polygons(list(Sr1), "s1")
Srs2 = Polygons(list(Sr2), "s2")
Srs3 = Polygons(list(Sr3), "s3")
Srs4 = Polygons(list(Sr4), "s4")

SpP = SpatialPolygons(list(Srs1,Srs2,Srs3,Srs4), 1:4)
plot(SpP, col = 1:4, pbg="white")
types <- letters[1:4]
polys <- SpatialPolygonsDataFrame(SpP,
                               data = data.frame(
                                 types=types,
                                 row.names = row.names(SpP)
                               ))
                               #extract


plot(r)
plot(polys,add=T)

functions_summ <-
function(x, na.rm = T){
c(
  mean = mean(x, na.rm = na.rm),
  min = min(x, na.rm = na.rm),
  max = max(x, na.rm = na.rm),
  sd = sd(x, na.rm = na.rm)
)
}

ex <- extract(stacked,polys,fun=functions_summ,df=T)
rownames(ex) <- 
paste0(gsub("\\..","",rownames(ex)),"_",rep(unique(polys@data$types),each=4))
#convert rownames to column
ex$var <- rownames(ex)

但我认为不可能将此data.frame附加到多边形的@data槽中,因为您拥有的观察值比多边形特征更多。您必须将数据调整为宽格式才能合并这两个数据库

票数 0
EN

Stack Overflow用户

发布于 2021-04-24 08:46:01

如果"Type“是多边形的属性,则这些多边形表示为行。使用单层,您可以

代码语言:javascript
运行
复制
library(terra)
v <- shapefile(system.file("external/lux.shp", package="raster"))
r <- raster(v, res=.1)
values(r) <- 1:ncell(r)

f <- function(x, na.rm = T) {
       c(mean=mean(x, na.rm = na.rm),
        min=min(x, na.rm = na.rm),
        max=max(x, na.rm = na.rm),
        sd=sd(x, na.rm = na.rm)
    )
}
e <- extract(r, v, f)

x <- cbind(v, e)

但这不适用于多个层(在RasterStack中);因此,您必须遍历这些层。

使用terra版本1.2-1 (当前为development version),您可以执行以下操作

代码语言:javascript
运行
复制
library(terra)
#terra version 1.2.1
r <- rast(system.file("ex/elev.tif", package="terra"))
x <- c(r, r/2)
names(x) <- c("a", "b")
v <- vect(system.file("ex/lux.shp", package="terra"))

f <- function(x, na.rm = T) {
       c(mean=mean(x, na.rm = na.rm),
        min=min(x, na.rm = na.rm),
        max=max(x, na.rm = na.rm),
        sd=sd(x, na.rm = na.rm)
    )
}

e <- extract(x, v, f)
head(e)
#  ID   a.mean a.min a.max     a.sd   b.mean b.min b.max     b.sd
#1  1 467.1052   339   547 34.58480 233.5526 169.5 273.5 17.29240
#2  2 333.8629   195   514 68.03058 166.9315  97.5 257.0 34.01529
#3  3 377.3712   256   517 77.14169 188.6856 128.0 258.5 38.57085
#4  4 373.6000   213   520 82.79129 186.8000 106.5 260.0 41.39564
#5  5 418.6490   293   511 48.35488 209.3245 146.5 255.5 24.17744
#6  6 314.9969   164   403 49.02558 157.4985  82.0 201.5 24.51279
票数 0
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/66205192

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档