[R-sig-Geo] Isolating only land areas on a global map for computing averages

[Tom Philippi]

The easiest approach would be to create a separate aligned raster layer for
land vs water.  There are plenty of coastline polygons available out there
(e.g., maptools, rworldmap, rnaturalearth packages): choose one in your
raster CRS (or choose one and spTransform() it).  Then, use a grid version
of your raster to extract values from that land/water SpatialPolygons
object.

1: Your idea of extracting the land/water value at each grid cell centroid
gives one estimate.  Instead of TRUE/FALSE, think of the numeric
equivalents 1,0,  then using those as weights for averaging across your
grid cells.
2: A "better" estimate would be to compute the fraction of each grid cell
that is land, then use those fractional [0, 1] values as weights to compute
a weighted average of precipitation over land.  At 2.8deg grid cells, a lot
of heavy rainfall coastal areas would have the grid cell centroid offshore
and be omitted by approach #1.
3: I recommend that you think hard about averaging across cells in Lat Lon
to estimate average precipitation over land.  The actual area of a ~2.8 by
2.8 deg grid cell at the equator is much larger than the same at 70 deg N.
I would spend the extra hour computing the actual area (in km^2) of land in
each of your 8192 grid cells, then using those in a raster as weights for
whatever calculations you do on the raster stack.  [Or you can cheat, as
the area of a grid cell in degrees is a function of only the latitudes, and
your required weights are multiplicative.]

Your mileage may vary...

Tom

On Wed, Nov 6, 2019 at 6:18 PM rain1290--- via R-sig-Geo <
r-sig-geo using r-project.org> wrote:

> Hi there,
> I am interested in calculating precipitation medians globally. However, I
> only want to isolate land areas to compute the median. I already first
> created a raster stack, called "RCP1pctCO2median", which contains the
> median values. There are 138 layers, with each layer representing one
> year.  This raster stack has the following attributes:
> class       : RasterStack
> dimensions  : 64, 128, 8192, 138  (nrow, ncol, ncell, nlayers)
> resolution  : 2.8125, 2.789327  (x, y)
> extent      : -181.4062, 178.5938, -89.25846, 89.25846  (xmin, xmax, ymin,
> ymax)
> coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
> names       :    layer.1,    layer.2,    layer.3,    layer.4,
> layer.5,    layer.6,    layer.7,    layer.8,    layer.9,   layer.10,
> layer.11,   layer.12,   layer.13,   layer.14,   layer.15, ...
> min values  : 0.42964514, 0.43375653, 0.51749371, 0.50838983, 0.45366730,
> 0.53099146, 0.49757186, 0.45697752, 0.41382199, 0.46082401, 0.45516687,
> 0.51857087, 0.41005131, 0.45956529, 0.47497867, ...
> max values  :   96.30350,  104.08584,   88.92751,   97.49373,
> 89.57201,   90.58570,   86.67651,   88.33519,   96.94720,  101.58247,
> 96.07792,   93.21948,   99.59785,   94.26218,   90.62138, ...
>
> Previously, I was isolating a specific region by specifying a range of
> longitudes and latitudes to obtain the medians for that region, like this:
> expansion1<-expand.grid(103:120, 3:15)lonlataaa<-extract(RCP1pctCO2Median,
> expansion1)Columnaaa<-colMeans(lonlataaa)
>
> However, with this approach, too much water can mix with land areas, and
> if I narrow the latitude/longitude range on land, I might miss too much
> land to compute the median meaningfully.
> Therefore, with this data, would it be possible to use an IF/ELSE
> statement to tell R that if the "center point" of each grid cell happens to
> fall on land, then it would be considered as land (i.e. that would be TRUE
> - if not, then FALSE)? Even if a grid cell happens to have water mixed with
> land, but the center point of the grid is on land, that would be considered
> land. But can R even tell what is land or water in this case?
> Thank you, and I would greatly appreciate any assistance!
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo using r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

	[[alternative HTML version deleted]]