2

I have two raster layers with the same projection, extent and resolution:

r1:

class       : RasterLayer 
dimensions  : 33434, 89431, 2990036054  (nrow, ncol, ncell)
resolution  : 0.001447053, 0.001447053  (x, y)
extent      : -56.50514, 72.90626, 24.28364, 72.66441  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs 
data source : 
names       : r1
values      : 11, 52  (min, max)
attributes  : ID    COUNT
        from: 11  9463473
        to  : 52 84036891

r2:

class       : RasterLayer 
dimensions  : 33434, 89431, 2990036054  (nrow, ncol, ncell)
resolution  : 0.001447053, 0.001447053  (x, y)
extent      : -56.50514, 72.90626, 24.28364, 72.66441  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs 
data source : C:\Users\arandaba\Desktop\test_gisR\webos.tif 
names       : r2
values      : 11, 99  (min, max)

Then I have a SpatialPolygonsDataFrame with same projection as above but different extent:

p:

class       : SpatialPolygonsDataFrame 
features    : 829 
extent      : -10.01233, 30.84247, 36.94062, 65.95678  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs 
variables   : 8
names       : SITE_CODE, ISO3, IUCNCAT, Year,  ss90,  sl90,  ss06,  sl06 
min values  :      1685,  AUT,       I, 1990,     0,     0,     0,     0 
max values  : 555552467,  SVN,      II, 2012, 65750, 82966, 68374, 89058 

When I extract raster values from those 829 polygons, I obtain a different number of observations for each raster layer. That is:

ext1 <- extract(r1,p,df=TRUE) # Results in 2040332 obs.
ext2 <- extract(r2,p,df=TRUE) # Results in 2040335 obs.

Since I want to examine changes in raster values from r1 to r2 within polygons, how could I solve this problem? Should I brick both raster layers and then run extract()? I need to do this efficiently, as I'm using large data files.

  • Welcome to GIS:SE sca! I edited the formatting of your question which looks quite different from before. If you don't agree then you can always re-edit your question and/or rollback to how you originally posted it :) – Joseph Apr 24 '15 at 10:23
  • Are you accounting for NA's? The default NA behavior for extract is na.rm=TRUE which would change the number of observations. – Jeffrey Evans Apr 24 '15 at 14:05
  • No, as I did not supply an argument fun. I just want to extract all pixel values within polygons. Thank you anyways! – sca Apr 24 '15 at 14:30
  • That is strange. What is your version of 'raster' (is it up to date?) extract(stack(r1, r2), p, df=TRUE) should get rid of this problem and should be more efficient. Creating a RasterBrick first would probably slow things down. – Robert Hijmans Apr 24 '15 at 16:23

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.