I am a novice in spatial statistics in R. So, I am looking for some advises. In particular, I would like to create a predictive map of capture success from a generalized linear mixed-effects model (GLMM). In my GLMM, the dependent variable is the total number of unique individuals captured per 100 trap nights in trapping sites and the independent variables are the proportion of land cover types in trapping sites. I included a random effect for the year. My raster of land cover types has cells of 25 x 25 m. How can I assign predicted values from the GLMM to raster pixels in R? Have I to use a kriging method? As the covariates are proportions of land cover types in trapping sites, is it a problem if a raster pixel represents a single land cover type (and not several land cover types)?


Here are:

  • The structure of my dataset:

enter image description here

  • A representation of one trapping site:

enter image description here

So, a trapping site is composed of several land cover types and each raster pixel is assigned to one land cover type. The trapping sites have different sizes. I defined the proportion of a given land cover type as number of pixels for the land cover type in the trapping site / total number of pixels (thus for all land cover types) in the trapping site.

  • Given that any well-defined trapping site is so small that it can be considered to have only one land cover type, exactly how did you define and measure "proportion of land cover types" at your sites? Did you look at the proportions within a certain distance, for instance?
    – whuber
    Feb 8, 2016 at 19:54
  • Thank you very much whuber for your response. No, I didn't look at the proportions within a certain distance, but within the entire trapping site.
    – Nell
    Feb 8, 2016 at 22:27
  • It looks like you have a change of support problem, especially if the trapping sites are of different sizes. In effect, your data represent properties of regions of relatively large spatial extent, whereas you wish to draw a map in which properties are assigned to points. Something has to give. If most of the trapping sites are approximately the same size, one work-around for this problem is to surround every grid pixel with a hypothetical trapping area of this common size (and shape) and assign the land cover proportions from the hypothetical area to the pixel.
    – whuber
    Feb 8, 2016 at 22:35
  • Thank you very much whuber. Yes, the trapping sites are of different sizes.
    – Nell
    Feb 8, 2016 at 22:47

1 Answer 1


If I understand your specified model correctly, a separate covariate for each landcover type proportion. Yes, the spatial data structure will need to be addressed. You will need to create a separate raster for each landcover type that represents proportion of each given type at the appropriate scale. Depending on how your independent variables were derived, this could take the form of a focal (NxN window or radius) or a zonal function (proportion within each discrete trapping unit).

I am going to mostly ignore any issues with model specification and aim you towards the answer. Although, I do encourage you to really think about your model structure and dig into fit before drawing inference or conclusion around a spatial estimate. Keep in mind that your independent variables are being drawn from the same process and, as such, are not entirely independent. This would be a case where something like AIC would really fail. I am also assuming that your random effect are the trapping units. Although, if you are controlling for the effect it would be a fixed effect. This is also most definitely a Modifiable Areal Unit Problem (MAUP) issue and you should dig into the literature to understand how this has effected your analysis. In your case, the most relevant literature will be in relation to analysis of game units (hunting zones).

The raster package provides an interface, for the raster object classes, to the generic predict function. The way the generic works is that, if needed, each package writes a predict function for the specific model object type produced by the package. This can be done in a way the when "predict" is called it recognizes the object class and calls the correct predict function. In this way you can just call predict on the raster object, representing the model covariates, and pass your model to it with the relevant arguments. This is a multilevel call because it first calls the raster predict and then that function calls the appropriate predict function for the model, all handled without you ever knowing (thanks Robert H).

Just keep in mind that your rasters need to represent your covariates. So, if your raster is nominal and not proportion of landcover types, then you will need to calculate a raster that represents the percents at the scale the model was constructed. If you just pulled the values of x from the raster than this is taken care of. Also, if your model was multivariate then your rasters need to share common number of rows, columns, origin coordinates and cell size. You would then read the data in as a raster stack or brick. I have a tutorial for implementing a random forests model using spatial data on my Quantitative Spatial Ecology website. This tutorial has example data and provides usage of the predict function for making a spatial estimate.

Now, here is where it may get ugly. I do not believe that the developers of glmm have written a formal predict function. Although, the lme4 library provides a predict function "predict.merMod" that should handle this object class. However, I am not sure that R's predict generic will catch that this is the predict function for this object class. Because of this you may need to explicitly pass the "fun" argument, in predict.raster, the "predict.merMod" funciton.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.