使用 R 创建分类栅格值的直方图？ （或者，使用 lat/long 值创建数据 table）
Use R to create a historgram of categorical raster values? (or, create data table with lat/long values)
总体而言，我对 R 和编程还很陌生，请原谅我即将到来的无能。
我正在处理大型分类栅格。本质上，大堡礁上每一个浅于 10 米的像素都被分配了一个值：11、12、13 或 15。（原始文件
here). My end goal is to create a histogram showing the frequency of the "rubble" category (which is the given by the value 12) by latitude. It would look very similar to the third panel of this figure，但他们显示“珊瑚栖息地”的地方我会显示瓦砾。
我认为最好的方法是尝试将原始栅格转换为数据框，其中每一行代表一个像素，并且有三列：分类值（11、12、13 或15)、纬度和经度。然后我可以使用这个数据框来创建任意数量的基本图。
理想情况下，我想在创建此数据框的过程中省略 NA，因为栅格为 152,505 x 112,421，但超过 99% 的像素是空的（由于大堡礁的形状）。
我可以轻松读取栅格并使用 Raster 或 Terra 包绘制它：
benthic <- ("data/GBR10 GBRMP Benthic.tif")
benthic_r <-rast(benthic)
plot(benthic_r)
我尝试使用它来缩小它，以便将来的计算更容易。
benthic2 <- na.omit(benthic)
benthic2_r <- rast(benthic2)
但是发现na.omit只适用于矢量数据，所以没有成功。
我正在尝试使用 Geocomputation with R 并想也许我会以某种方式使用 as.data.frame 函数，但我一直无法找到一种方法来使用它来创建 table我要。
我也一直在浏览有关 Stars 和 RasterVis 包的信息。我试图跳过数据框的想法并使用 rasterVis 包直接从原始栅格创建直方图：
dirY <- xyLayer(benthic_r, y) #trying to assign y axis to be the y values of the pixels (which are hopefully latitude...)
dirXY <- xyLayer(count(12), benthic_r) #trying to assign x values to be a count of the number of pixels with a value of 12
histogram(benthic_r, dirXY,
maxpixels = 1e+05,
strip=TRUE,
par.settings=rasterTheme(),
att = 1)
#attempting to follow https://cran.r-project.org/web/packages/rasterVis/rasterVis.pdf "histogram methods" section
我在很多方面对此进行了调整，但一切都会导致新的错误。
我不期待完整的代码解决方案，而是如果您能帮助我了解在哪里查看以及如何逐步构建以创建该图表，我将不胜感激。首先尝试创建数据框是一种好方法吗？还是完全没有必要？
我找到的关于如何处理栅格数据的资源似乎都是关于处理卫星图像（连续栅格）的，我没有找到很多处理分类栅格的技术。如果您知道有关处理分类栅格的优秀教程，那将不胜感激。
提前感谢您的任何建议。
我认为最简单的方法是聚合栅格，以便获得一列，其中包含感兴趣的值。也就是说，每个纬度（行）对应一个单元格。这是一个例子：
# example data
library(terra)
r <- rast(nrow=18)
values(r) <- sample(c(11,12,13,15), ncell(r), replace=TRUE)
# count the number of cells that are 12
r12 <- r == 12
# aggregate the columns
a <- aggregate(r12, c(1, ncol(r)), sum)
# get the latitude for each row
lat <- yFromRow(a, 1:nrow(a))
您可以制作不同类型的图。这里有两个例子。我更喜欢纵轴上的纬度。
plot(values(a), lat, ylab="latitude", xlab="count")
barplot(rev(values(a)[,1]), horiz=T, names.arg=rev(lat), las=1)
总体而言，我对 R 和编程还很陌生，请原谅我即将到来的无能。
我正在处理大型分类栅格。本质上，大堡礁上每一个浅于 10 米的像素都被分配了一个值：11、12、13 或 15。（原始文件 here). My end goal is to create a histogram showing the frequency of the "rubble" category (which is the given by the value 12) by latitude. It would look very similar to the third panel of this figure，但他们显示“珊瑚栖息地”的地方我会显示瓦砾。
我认为最好的方法是尝试将原始栅格转换为数据框，其中每一行代表一个像素，并且有三列：分类值（11、12、13 或15)、纬度和经度。然后我可以使用这个数据框来创建任意数量的基本图。
理想情况下，我想在创建此数据框的过程中省略 NA，因为栅格为 152,505 x 112,421，但超过 99% 的像素是空的（由于大堡礁的形状）。
我可以轻松读取栅格并使用 Raster 或 Terra 包绘制它：
benthic <- ("data/GBR10 GBRMP Benthic.tif")
benthic_r <-rast(benthic)
plot(benthic_r)
我尝试使用它来缩小它，以便将来的计算更容易。
benthic2 <- na.omit(benthic)
benthic2_r <- rast(benthic2)
但是发现na.omit只适用于矢量数据，所以没有成功。
我正在尝试使用 Geocomputation with R 并想也许我会以某种方式使用 as.data.frame 函数，但我一直无法找到一种方法来使用它来创建 table我要。
我也一直在浏览有关 Stars 和 RasterVis 包的信息。我试图跳过数据框的想法并使用 rasterVis 包直接从原始栅格创建直方图：
dirY <- xyLayer(benthic_r, y) #trying to assign y axis to be the y values of the pixels (which are hopefully latitude...)
dirXY <- xyLayer(count(12), benthic_r) #trying to assign x values to be a count of the number of pixels with a value of 12
histogram(benthic_r, dirXY,
maxpixels = 1e+05,
strip=TRUE,
par.settings=rasterTheme(),
att = 1)
#attempting to follow https://cran.r-project.org/web/packages/rasterVis/rasterVis.pdf "histogram methods" section
我在很多方面对此进行了调整，但一切都会导致新的错误。
我不期待完整的代码解决方案，而是如果您能帮助我了解在哪里查看以及如何逐步构建以创建该图表，我将不胜感激。首先尝试创建数据框是一种好方法吗？还是完全没有必要？
我找到的关于如何处理栅格数据的资源似乎都是关于处理卫星图像（连续栅格）的，我没有找到很多处理分类栅格的技术。如果您知道有关处理分类栅格的优秀教程，那将不胜感激。
提前感谢您的任何建议。
我认为最简单的方法是聚合栅格，以便获得一列，其中包含感兴趣的值。也就是说，每个纬度（行）对应一个单元格。这是一个例子：
# example data
library(terra)
r <- rast(nrow=18)
values(r) <- sample(c(11,12,13,15), ncell(r), replace=TRUE)
# count the number of cells that are 12
r12 <- r == 12
# aggregate the columns
a <- aggregate(r12, c(1, ncol(r)), sum)
# get the latitude for each row
lat <- yFromRow(a, 1:nrow(a))
您可以制作不同类型的图。这里有两个例子。我更喜欢纵轴上的纬度。
plot(values(a), lat, ylab="latitude", xlab="count")
barplot(rev(values(a)[,1]), horiz=T, names.arg=rev(lat), las=1)