rareawgs84r-rastermaxent

Calculating area of occupancy from a binary unprojected raster


I have a series of binary raster layers (ascii file) showing presence/absence of a species in Europe and Africa. The file is based on unprojected lat/long (WGS84) data. My aim is to calculate the area of presence using R (I don't have access to ArcGIS).

I know that the raster package has a function for calculating area, but I'm worried that this won't be accurate for unprojected data. I have also looked at the cellStats function in the raster package, and can use this to "sum" the number of cells occupied, but I feel this has the same problem.

jan<-raster("/filelocation/file.asc")
jan
class       : RasterLayer 
dimensions  : 13800, 9600, 132480000  (nrow, ncol, ncell)
resolution  : 0.008333333, 0.008333333  (x, y)
extent      : -20, 60, -40, 75  (xmin, xmax, ymin, ymax)
coord. ref. : NA 
data source : "/filelocation"
names       : file.asc
values      : -2147483648, 2147483647  (min, max)

area(jan)
class       : RasterLayer 
dimensions  : 13800, 9600, 132480000  (nrow, ncol, ncell)
resolution  : 0.008333333, 0.008333333  (x, y)
extent      : -20, 60, -40, 75  (xmin, xmax, ymin, ymax)
coord. ref. : NA 
names       : layer 
values      : 6.944444e-05, 6.944444e-05  (min, max)

Warning messages:
1: In .local(x, ...) :
  This function is only useful for Raster* objects with a longitude/latitude     coordinates
2: In .rasterFromRasterFile(grdfile, band = band, objecttype, ...) :
  size of values file does not match the number of cells (given the data type)

cellStats(jan,"sum")
[1] 3559779

Anybody know of a way to calculate the presence area accurately, accounting for the earth curvature?

Thanks!


Solution

  • I do not know what is going in with your file (why you get warning #2). But is here is a work around

    r <- raster(nrow=13800, ncol=9600, xmn=-20, xmx=60, ymn=-40, ymx=75)
    # equivalent to r <- raster(jan)
    x = area(r)
    x
    class       : RasterLayer 
    dimensions  : 13800, 9600, 132480000  (nrow, ncol, ncell)
    resolution  : 0.008333333, 0.008333333  (x, y)
    extent      : -20, 60, -40, 75  (xmin, xmax, ymin, ymax)
    coord. ref. : +proj=longlat +datum=WGS84 
    data source : c:\temp\R_raster_Robert\2015-01-26_213612_1208_85354.grd 
    names       : layer 
    values      : 0.2227891, 0.8605576  (min, max)
    

    Now you have the area of each cell in km2. By multiplying these values with Raster objects with presence/absence values and then using cellStats( , 'sum') you can obtain the total area with presence.