I am very new at using GIS data and only modestly experienced with R. I've been reading about how to analyze spatial data using the spatial-analyst.net PDF book, so I'm not completely lost, but I thought I could describe my problem and people might suggest ideas.

I have a dataset with about 2000 measurements at different lat/long coordinates, although I will probably subdivide this dataset as the data was collected over 3 years and conditions changed over time. Let's call the variable being measured "IP."

I want to create a map of IP in the full area under question using Kriging or some other interpolation method on the sample data. Then I want to create a histogram measuring the amount of land in various IP buckets. I will also need to create a histogram which shows the number of samples in each bucket (note a sample could have a higher or lower actual IP than what kriging predicts for its land).

I follow how to load the data into a SpatialPointsDataFrame and run a kriging analysis, where I'm having trouble is how to convert that data into a gridded dataframe so I can do the histogram analysis.

Any suggestions for converting points into grids?

3 Answers 3


You're right ... it is pretty easy! The "raster" package has some pretty straightforward ways of dealing with creating and manipulating rasters.


# Load your point shapefile (with IP values in an IP field):
pts <- readShapePoints("pts.shp")

# Create a raster, give it the same extent as the points
# and define rows and columns:

rast <- raster()
extent(rast) <- extent(pts) # this might be unnecessary
ncol(rast) <- 20 # this is one way of assigning cell size / resolution
nrow(rast) <- 20

# And then ... rasterize it! This creates a grid version 
# of your points using the cells of rast, values from the IP field:
rast2 <- rasterize(pts, rast, pts$IP, fun=mean) 

You can assign grid size and resolution in a number of ways - have a good look at the raster package documentation.

The values of the raster cells from rasterize can be calculated with a function - 'mean' in the example above. Make sure you put this in: otherwise it just uses the value of IP from the last point it comes across!

From a CSV:

pts <- read.csv("IP.csv")
coordinates(pts) <- ~lon+lat
rast <- raster(ncol = 10, nrow = 10)
extent(rast) <- extent(pts)
rasterize(pts, rast, pts$IP, fun = mean)
  • Hey this is very useful, but how would the code look if I started with the points in a simple CSV with lat/longs rather than a shapefile? The columns in the CSV would be IP, Lat, Long, etc, etc, etc. May 2, 2012 at 17:03
  • You did indicate that you've already loaded the data into a SpatialPointsDataFrame ... which is exactly what pts is in my example above. Just run the code on your SpatialPointsDataFrame object!
    – Simbamangu
    May 2, 2012 at 18:33
  • 4
    This reply, although excellent, does not seem to address what is needed. (It appears to offer a solution to gis.stackexchange.com/questions/20018 instead.) The challenge is to interpolate 2000 or so points, not just assign their values to raster cells. Given that the OP claims already to have "run a kriging analysis," this question comes down to extracting the values of a raster (say, r) for using in a hist-like procedure, which is simply a matter of an expression like hist(getValues(r)).
    – whuber
    May 2, 2012 at 19:46
  • @whuber - Looks like OP asks "where I'm having trouble is how to convert that data into a gridded dataframe so I can do the histogram analysis ... any suggestions for converting points into grids" as the actual question, and knows how to make a SpatialPointsDataFrame and run the kriging. But you're right, it seems to be a duplicate of 20018 (except for gridded input).
    – Simbamangu
    May 3, 2012 at 3:37
  • Apologies, @user1080253 ... I read 'grid' as 'raster' which isn't correct, and not helpful for kriging; see here for a better idea on creating a regular grid and interpolating your data to that grid.
    – Simbamangu
    May 3, 2012 at 3:51

The plotKML package has a function called vect2rast. This function basically extends the rasterize function available in the raster package. The advantage of vect2rast; however, is that it requires no input from the user's side i.e. it automatically determines the grid cell size and the bounding box based on the properties of the input data set. The grid cell size is estimated based on the density/size of features in the map (nndist function in spatstat package).

Rast2 <- vect2rast(pts)

# for large data sets use SAGA GIS:
Rast2 <- vect2rast(pts, method = "SAGA")

If you are starting from a spatialPoints object you just do the following:

# Convert spatialpoints object to a 'conventional' dataframe
spatial_points_dataframe  <- data.frame(spatial_points_object)

# View the first few rows

Your Answer

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

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