Search code examples
rmapslevelplot

Add XY points conditionally to raster map generated by levelplot


I have a tricky ifelse task here. Below is my code showing data for two periods(future and current). The data has mean, 5th and 95th confidence bounds alongside the xy cordinates. I want to compare the mean, 5th and 95th confidence bounds (CI) for the two dfs (future and current).

Conditions:

1) if the CIs for future do not overlap with those of current and CIs of futrue are above current, thenpch=2.

2) if the CIs for future do not overlap with those of current and CIs of futrue are below current, thenpch=3.

3) if CIs of future overlap with those of current, then pch=4

library(raster)
library(rasterVis)

    s <- stack(replicate(2, raster(matrix(runif(100), 3))))
    current <- data.frame(coordinates(sampleRandom(s, 3, sp=TRUE)),
                     C5th=c(17.643981,16.83572,9.979904),
                     CMean=c(26.66364,19.74286,15.10000),C95th=c(35.68329,22.64999,20.22010))


    future <- data.frame(coordinates(sampleRandom(s, 3, sp=TRUE)),
                          C5th=c(17.643981,16.83572,9.979904)*2,
                          CMean=c(26.66364,19.74286,15.10000)*2,C95th=c(35.68329,22.64999,20.22010)*2)

The result of the above three conditions will then be added to my map. Something like (just an attempt):

levelplot(s, margin=FALSE, at=seq(0, 1, 0.05)) + 
  layer(sp.points(xy, pch=ifelse(condition, 2, 3,4), cex=2, col=1), columns=1) +
  layer(sp.points(xy, pch=ifelse(condition, 2, 3,4), cex=2, col=1), columns=2)

For example, in the figure below, if the minimum of NFC (future) lies completely above the maximum of AFC (current), then condtion 1. if the maximum of NFC lies completely below the minimum of AFC, then condtion 1.The plot as shown below satisfies condition 3.

enter image description here

Please help. AT.


Solution

  • It is easier if you define a complete SpatialPointsDataFrame object with an additional categorical variable defined according to the conditions you need.

    library(raster)
    library(rasterVis)
    
    s <- stack(replicate(2, raster(matrix(runif(1000), 3))))
    ## Coordinates
    cc <- sampleRandom(s, 3, sp = TRUE)
    ## Data
    current <- data.frame(C5th=c(17.643981,16.83572,9.979904),
                          CMean=c(26.66364,19.74286,15.10000),
                          C95th=c(35.68329,22.64999,20.22010))
    
    future <- data.frame(C5th=c(17.643981,16.83572,9.979904)*2,
                         CMean=c(26.66364,19.74286,15.10000)*2,
                         C95th=c(35.68329,22.64999,20.22010)*2)
    
    cf <- data.frame(current, future)
    ## Define a categorical variable checking the conditions you need
    cf$class <- with(cf,
                     ifelse(C5th > C95th.1, 'A',
                            ifelse(C95th < C5th.1, 'B',
                                   ifelse(C5th < C95th.1 && C5th > C5th.1, 'C', 'D')
                                   )
                            )
                     )
    cf$class <- factor(cf$class)
    
    ## Build a SPDF object with the coordinates and data
    pp <- SpatialPointsDataFrame(cc, cf)
    

    This object can be displayed with spplot. With it you can choose the symbols, sizes, etc.

    levelplot(s) + spplot(pp["class"],
                          pch = 21:23,
                          col.regions = 'gray')