r почему извлечение растра дает тусклую ошибку?

Я пытаюсь извлечь значения и получить среднее значение из растрового блока, но получаю ошибку, которая, как мне кажется, связана с размерами растрового блока.

данные были загружены с сайта NOAA

я сделал следующее:

library(raster)

ERSST <- rotate(brick('sst.mnmean.nc'))
ERSST
class       : RasterBrick 
dimensions  : 89, 180, 16020, 1976  (nrow, ncol, ncell, nlayers)
resolution  : 2, 2  (x, y)
extent      : -179, 181, -89, 89  (xmin, xmax, ymin, ymax) # ignore extent, needs fixing but not relevant for the question
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
data source : in memory
names       : X1854.01.01, X1854.02.01, X1854.03.01, X1854.04.01, X1854.05.01, X1854.06.01, X1854.07.01, X1854.08.01, X1854.09.01, X1854.10.01, X1854.11.01, X1854.12.01, X1855.01.01, X1855.02.01, X1855.03.01, ... 
min values  :        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8,        -1.8, ... 
max values  :    32.09937,    31.44189,    31.72137,    31.47466,    33.23633,    31.90788,    35.22922,    33.66898,    32.26702,    32.15502,    31.68270,    31.74512,    31.32458,    31.23049,    29.88974, ... 
Date        : 1854-01-01, 2018-08-01 (min, max)

Произвольная точка

xy <- data.frame(x = -49, y = 45)

и когда я извлекаю, я получаю:

extract(ERSST, xy, buffer = 1e+05, small = TRUE, fun = mean)
Error in apply(x, 2, fun2) : dim(X) must have a positive length

что заставляет меня думать, что это проблема измерения, так это ошибка, которую я получаю, когда пытаюсь указать слои для использования

extract(ERSST, xy, buffer = 1e+05, small = TRUE, layer = 10, nl = 10)
Error in x[, lyrs] : incorrect number of dimensions

и, кажется, работает нормально, если я сначала усредняю (но это не то, что я хочу, мне нужен временной ряд в данный момент)

mERSST <- mean(ERSST)
extract(mERSST, xy, buffer = 1e+05, small = TRUE, fun = mean)
[1] 5.649212

Возможно, это атрибут Date в растровом блоке. Любые обходные пути или решения для предотвращения этой ошибки?

Ответ @RobertHijmans заставил меня понять, что я всегда получаю только одно значение из извлечения, даже когда точка находится на стыке нескольких ячеек сетки, как в примере выше.

plot(MSST, xlim = c(-60, -40), ylim = c(40, 50))
points(xy)

С использованием:

extract(mERSST, xy, buffer = 1e+05, small = TRUE, cellnumbers = TRUE)
[[1]]
       cell       value 
3845.000000    5.649212

Я получаю только одно значение, тогда как я ожидал бы, что их будет 4, независимо от того, насколько мал буфер. Я что-то пропустил в выписке? Поэтому я попытался преобразовать свою точку в круг и использовать ее для извлечения данных.

coordinates(xy) <- ~ x + y
proj4string(xy) <- '+init=epsg:4326'

xy_utm <- spTransform(xy, CRS('+init=epsg:32621'))
gbf_utm <- rgeos::gBuffer(xy_utm, width = 1e5, quadsegs = 250L)
gbf <- spTransform(gbf_utm, CRS(proj4string(xy)))


plot(ERSST[[1]], xlim = c(-60, -40), ylim = c(40, 50))
points(xy, pch = 19)
plot(gbf, add = TRUE)

extract(ERSST[[1]], gbf, small = TRUE, weights = TRUE)

это дает мне:

[[1]]
        value weight
[1,] 1.722664   0.25
[2,] 3.683457   0.25
[3,] 5.985203   0.25
[4,] 8.442450   0.25

в версии 2.6.7 (и это, кажется, имеет смысл).

но

[[1]]
          value      weight
  [1,] 1.722664 0.001236928
  [2,] 1.722664 0.003935680
  [3,] 1.722664 0.005285056
  [4,] 3.683457 0.005285056
  [5,] 3.683457 0.003935680
  [6,] 3.683457 0.001236928
  [7,] 1.722664 0.002136512
  [8,] 1.722664 0.008321151
  [9,] 1.722664 0.011244799
 [10,] 1.722664 0.011244799
 [11,] 1.722664 0.011244799
 [12,] 3.683457 0.011244799
 [13,] 3.683457 0.011244799
 [14,] 3.683457 0.011244799
 [15,] 3.683457 0.008208703
 [16,] 3.683457 0.001911616
 [17,] 1.722664 0.003036096
 [18,] 1.722664 0.010907455
 [19,] 1.722664 0.011244799
 [20,] 1.722664 0.011244799
 [21,] 1.722664 0.011244799
 [22,] 1.722664 0.011244799
 [23,] 3.683457 0.011244799
 [24,] 3.683457 0.011244799
 [25,] 3.683457 0.011244799
 [26,] 3.683457 0.011244799
 [27,] 3.683457 0.010907455
 [28,] 3.683457 0.003036096
 [29,] 1.722664 0.000449792
 [30,] 1.722664 0.010232767
 [31,] 1.722664 0.011244799
 [32,] 1.722664 0.011244799
 [33,] 1.722664 0.011244799
 [34,] 1.722664 0.011244799
 [35,] 1.722664 0.011244799
 [36,] 3.683457 0.011244799
 [37,] 3.683457 0.011244799
 [38,] 3.683457 0.011244799
 [39,] 3.683457 0.011244799
 [40,] 3.683457 0.011244799
 [41,] 3.683457 0.010232767
 [42,] 3.683457 0.000337344
 [43,] 1.722664 0.003036096
 [44,] 1.722664 0.011244799
 [45,] 1.722664 0.011244799
 [46,] 1.722664 0.011244799
 [47,] 1.722664 0.011244799
 [48,] 1.722664 0.011244799
 [49,] 1.722664 0.011244799
 [50,] 3.683457 0.011244799
 [51,] 3.683457 0.011244799
 [52,] 3.683457 0.011244799
 [53,] 3.683457 0.011244799
 [54,] 3.683457 0.011244799
 [55,] 3.683457 0.011244799
 [56,] 3.683457 0.002923648
 [57,] 5.985203 0.002923648
 [58,] 5.985203 0.011244799
 [59,] 5.985203 0.011244799
 [60,] 5.985203 0.011244799
 [61,] 5.985203 0.011244799
 [62,] 5.985203 0.011244799
 [63,] 5.985203 0.011244799
 [64,] 8.442450 0.011244799
 [65,] 8.442450 0.011244799
 [66,] 8.442450 0.011244799
 [67,] 8.442450 0.011244799
 [68,] 8.442450 0.011244799
 [69,] 8.442450 0.011244799
 [70,] 8.442450 0.002923648
 [71,] 5.985203 0.000337344
 [72,] 5.985203 0.010120319
 [73,] 5.985203 0.011244799
 [74,] 5.985203 0.011244799
 [75,] 5.985203 0.011244799
 [76,] 5.985203 0.011244799
 [77,] 5.985203 0.011244799
 [78,] 8.442450 0.011244799
 [79,] 8.442450 0.011244799
 [80,] 8.442450 0.011244799
 [81,] 8.442450 0.011244799
 [82,] 8.442450 0.011244799
 [83,] 8.442450 0.010007871
 [84,] 8.442450 0.000224896
 [85,] 5.985203 0.002811200
 [86,] 5.985203 0.010795007
 [87,] 5.985203 0.011244799
 [88,] 5.985203 0.011244799
 [89,] 5.985203 0.011244799
 [90,] 5.985203 0.011244799
 [91,] 8.442450 0.011244799
 [92,] 8.442450 0.011244799
 [93,] 8.442450 0.011244799
 [94,] 8.442450 0.011244799
 [95,] 8.442450 0.010682559
 [96,] 8.442450 0.002698752
 [97,] 5.985203 0.001799168
 [98,] 5.985203 0.007871359
 [99,] 5.985203 0.011244799
[100,] 5.985203 0.011244799
[101,] 5.985203 0.011244799
[102,] 8.442450 0.011244799
[103,] 8.442450 0.011244799
[104,] 8.442450 0.011244799
[105,] 8.442450 0.007871359
[106,] 8.442450 0.001799168
[107,] 5.985203 0.001236928
[108,] 5.985203 0.003935680
[109,] 5.985203 0.005285056
[110,] 8.442450 0.005285056
[111,] 8.442450 0.003935680
[112,] 8.442450 0.001236928

в версии 2.7-13, что не может быть правильным.

r raster r-raster

Lukas 14.09.2018 источник

comment

Это ошибка. Бывает, когда буфер мал такой, что выделена только одна ячейка (и внутренний результат — вектор, а не матрица, такая, что применить уже не получится). Я исправил это в версии: 2.7-11. Он доступен здесь: github.com/rspatial/raster - Robert Hijmans 15.09.2018

comment

@RobertHijmans, спасибо. Однако меня удивляет, что извлечение всегда возвращает одно значение, даже когда я пытаюсь извлечь значение прямо на границе двух (или четырех) ячеек сетки. Это почему? Я ожидаю значение 1, даже если буфер мал. - Lukas 17.09.2018

comment

При отсутствии буфера точка может находиться только в одной ячейке. Если пойдет направо и вниз, если на границе. Жизнь станет очень сложной, если вы не можете заранее знать, сколько значений будет возвращено. - Robert Hijmans 18.09.2018

comment

@RobertHijmans Достаточно справедливо для точки без буфера, но в приведенном выше примере есть буфер в 100 км, который должен покрывать ~ 4 ячейки. - Lukas 18.09.2018

comment

@RobertHijmans: я пытался указать многоугольник. Кажется, это работает в версии 2.6.7, где я получаю 4 значения и веса, которые имеют смысл, но когда я делаю то же самое в версиях 2.7-11 или 2.7-13, я получаю 112 значений. - Lukas 18.09.2018

comment

В вашей первой точке --- дополнительные ячейки включаются только в том случае, если буфер достигает центра ячейки. - Robert Hijmans 18.09.2018

comment

По вашему второму пункту. Это некрасиво, и я исправлю это, спасибо за отчет --- но обратите внимание, что веса по-прежнему верны, если вы их суммируете tapply(x$weight, x$value, sum) - Robert Hijmans 18.09.2018

Ответы (2)

arrow_upward
1
arrow_downward

Чтобы проиллюстрировать приведенное выше обсуждение в упрощенной форме (и после исправления проблемы с извлечением с полигонами), я получаю

library(raster)

r <- raster(xmn=-59, xmx=-39, ymn=41, ymx=49, res=2, vals=1:40)

xy <- SpatialPoints(data.frame(x = -49, y = 45), proj4string = CRS('+init=epsg:4326'))
p <- buffer(xy, width = 1e5, quadsegs = 250L)
plot(r)
plot(p, add=T)

extract(r, xy)       
#26 
extract(r, p)
#[[1]]
#[1] 16 26 25 15

extract(r, p, weights=T)
#[[1]]
#     value weight
#[1,]    15   0.25
#[2,]    16   0.25
#[3,]    25   0.25
#[4,]    26   0.25

extract(r, xy, buffer=100000)
#[[1]]
#value 
#   15 

extract(r, xy, buffer=1000000)
#[[1]]
# [1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
#[37] 37 38 39 40

Robert Hijmans 18.09.2018

arrow_upward
1
arrow_downward

Похоже, аргументы layer и nl не работают. Я не уверен, почему.

Один из способов обхода — сначала извлечь значения, а затем подмножество значений.

library(raster)

value <- extract(ERSST, point, buffer = 1e+05, small = TRUE)
value[[1]][10:19]
# X1854.10.01 X1854.11.01 X1854.12.01 X1855.01.01 X1855.02.01 X1855.03.01 X1855.04.01 X1855.05.01 
# 7.932416    5.712043    4.428292    3.010927    2.289096    2.385752    2.528488    3.261783 
# X1855.06.01 X1855.07.01 
# 5.762860    7.617740

www 14.09.2018

r почему извлечение растра дает тусклую ошибку?

Ответы (2)

Похожие вопросы