Take the 2-minute tour ×

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?

share|improve this question

1 Answer 1

up vote 7 down vote accepted

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) <- ~lat+lon
rast <- raster(ncol = 10, nrow = 10)
extent(rast) <- extent(pts)
share|improve this answer
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. –  user1080253 May 2 '12 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 '12 at 18:33
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 '12 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 '12 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 '12 at 3:51

Your Answer


