[R-sig-Geo] Raster Layers same resolution but not the same coordinates over the same area

Justin Michell jwm302 at gmail.com
Wed Jun 25 20:54:47 CEST 2014


Dear Group

I am attempting to plot the length of malaria risk (in months) over Kenya (using the model described in the paper by Gemperli et al. (2003) entitled, 'Mapping malaria transmission in West and Central Africa’). The criteria for the model as it pertains to which months are suitable for stable malaria transmission at the locations I have in my sample, can be viewed here: https://db.tt/AoThPoSC

Following these criteria l have obtained which months are suitable at each of the locations in my dataset, as follows:

# evaluate month 1 (January)
> df$month1 <- 1
df[!(df$MATemp1 > (19.5 + df$SDMeanTemp)  & (df$meanRain1>80 || df$meanRain12 > 80 || df$meanRain11>80) &
(df$MARain1 >  60 || df$NDVI12  > 0.35) & df$minAnnualTemp > 5), "month1"] <- 0

My problem is that now I want to use the same logic to obtain a malaria seasonal map covering the whole country, not just evaluating whether a month is suitable at my locations but everywhere that I have raster data.

I want to output something like this but for Kenya: https://db.tt/x0izxff7

I downloaded my NDVI data from the  Land Processes Distributed Active Archive Center (LP DAAC) (https://lpdaac.usgs.gov/, product MOD13A3), and the rest from the WorldClim website. All layers are at 1 km resolution and the NDVI data was resampled to the same geographic projection as temp and rain layers (I used the raster and gdalutils packages to get the data in and for processing). 

My raster layers are in stacks, one stack per enviro/climate layer made up of 12 months each. So meanTempStack has 12 layers representing an aggregated raster layer for each month of the year. I have converted these stacks to data frames keeping track of the coordinates.  These stack data frames (which are in a .rda workspace file) can be downloaded from: 

https://www.dropbox.com/s/uudzadx811gvnv2/rasterStacks.rda

And then once the .rda file is in your workspace, they can be retrieved in the usual way using:

> load("rasterStacks.rda”)

My plan was to merge all these data frames and then apply the same logic as demonstrated above but now for everywhere in Kenya not just at my locations. The problem can be demonstrated by the following (x denotes longitude y latitude):

> merge2 <- merge(meanTempStackDf, rainStackDf,  by = c("x", "y”))
> head(merge2)
        x            y meanTemp1 meanTemp2 meanTemp3 meanTemp4 meanTemp5 meanTemp6 meanTemp7 meanTemp8 meanTemp9 meanTemp10 meanTemp11 meanTemp12 rain1 rain2 rain3 rain4 rain5 rain6 rain7 rain8 rain9 1 33.9125 -0.004166667      23.1      23.3      23.4      22.9      22.4      21.9      21.6      21.6      22.1       22.7       22.6       22.5    52    68   129   172   130    44    34    46    51 2 33.9125  -0.012500000      23.1      23.3      23.4      22.9      22.4      21.9      21.6      21.6      22.1       22.7       22.6       22.5    52    68   128   172   130    44    34    45    50 3 33.9125 -0.020833333      23.1      23.3      23.4      22.9      22.4      21.9      21.6      21.6      22.1       22.7       22.6       22.5    52    68   128   171   130    44    34    45    50 4 33.9125 -0.029166667      23.1      23.3      23.4      22.9      22.4      21.9      21.6      21.6      22.1       22.7       22.7       22.5    51    68   128   171   130    44    34    44    49 5 33.9125 -0.037500000      23.1      23.3      23.4      22.9      22.4      21.9      21.6      21.6      22.1       22.7       22.7       22.5    51    68   128   170   130    44    34    44    49
6 33.9125 -0.045833333      23.1      23.3      23.4      22.9      22.4      21.9      21.6      21.6      22.1       22.7       22.7       22.5    51    68   128   170   130    44    33    44    48  rain10 rain11 rain12 1     73    119     85 2     73    119     85 3     72    119     85 4     72    118     86 5     71    118     86 6     71    118     86

But when I try to merge NDVI with a WorlClim layer there is no ‘overlapping’ of coordinates:

> merge2 <- merge(meanTempStackDf, rainStackDf,  by = c("x", "y”))
> head(merge2)
 [1] x      y      NDVI1  NDVI2  NDVI3  NDVI4  NDVI5  NDVI6  NDVI7  NDVI8  NDVI9  NDVI10 NDVI11 NDVI12 rain1  rain2  rain3  rain4  rain5  rain6  rain7  rain8  rain9  rain10 rain11 rain12
<0 rows> (or 0-length row.names)

So I am not exactly sure what to do from here. Any thoughts/help would be appreciated. 

Thanks
Justin Michell


More information about the R-sig-Geo mailing list