raster - Converting point data into gridded dataframe for histogram analysis using R?

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?

Answer

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

library(maptools)
library(raster)

# 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)