Question re: Intamap package / Interpolate method: Inputs NOT both Spatial Points Data Frame?
Thanks for pointing out that packages intamap and raster make incompatible assumptions about the function (in raster: method) interpolate. Sth we need to sort out. In your script, you assign gddThisYear$value <- gddThisYear[,jCtr] it is too bad that this does not give a warning or error, apparently one can assign S4 objects to data.frame columns -- it might even be documented behaviour. anyway, check the class of gddThisYear[,jCtr] yourself and try instead gddThisYear$value <- gddThisYear[[jCtr]] Hth,
On 10/26/2010 10:43 PM, Rick Reeves wrote:
Hello Edzer: (by the way, I am using R 2.12) Yes, I had - BUT - I rechecked my script, and the 'interpolate' method in intamap was apparently MASKED by the 'interpolate' method in the subsequently loaded raster package. (in-progress code is below...) I canceled the loading of the raster package, and re-ran the script with better results. This time, the line: interpImage <- interpolate(gddThisYear,monarchStaLocs,list(mean=TRUE, variance=TRUE)) results in: ------------------------------------------------- R 2010-10-26 13:31:48 interpolating 133 observations, 94 prediction locations [1] "estimated time for copula 60.5284609039258" Error in min(dataObs) : invalid 'type' (S4) of argument In addition: Warning message: In predictTime(nObs = dim(observations)[1], nPred = dim(coordinates(predictionLocations))[1], : using standard model for estimating time. For better platform spesific predictions, please run timeModels <- generateTimeModels() and save the workspace ------------------------------------------------ So intamap::interpolate documentation matches the argument list. But, my script is still not correct. I suspect that this error happens because I do not fully understand how to use 'interpolate' for this task. My task could also be solved by an approach like that demonstrated in Dylan Beaudette's excellent example at: http://casoilresource.lawr.ucdavis.edu/drupal/node/442 .. but using intamap::interpolate would be much simpler if I have made the correct assumptions. What are your thoughts? Thanks, Rick R ############################################################################################################## interpGddData <- function() { # require(gstat) require(intamap) # require(rgdal) require(maptools) # require(raster) setwd("v:/MigrationDynamics/NCEAS PROCESS") firstGddCol <- 5 # setwd("/data/scientist_share/scientist/MigrationDynamics/NCEAS PROCESS") gddTable <- readShapePoints("OhioGDDStationsGP.shp") monarchStaLocs <- readShapePoints("OH_Locations.shp") numGddCols <- 1 + (ncol(gddTable) - firstGddCol) # The Idea: For each sampling year, generate one map for each GDD parameter. # Each map is an interpolated 'surface' from all of the GDD-calculated 'station' # locations in the input file. theYears <- unique(gddTable at data$YR) numYears <- length(theYears) for (curYear in theYears) { gddThisYear <- gddTable[which(gddTable at data$YR==curYear),] # table entries: all stations, current year # make input table copy to append with interpolated GDD values for current year monarchStasWithGDD <- gddThisYear # Write one output ESRI Shape File for each year, containing the Monarch station # location points with all of the interpolated GDD parameters for the current year appended. outFileName <- sprintf("OH_MonarchSta_GDD_%d.tif",curYear) # theBase <- [,1:(firstGddColumn-1)] # Create one interpolated image for each parameter, resulting in a named matrix file. for (jCtr in firstGddCol : numGddCols) { # fieldName <- names(gddThisYear)[jCtr] # we need this for the fileName # fileName <- sprintf("%s_%d.tif",fieldName,curYear) gddThisYear$value <- gddThisYear[,jCtr] # Note: gddThisYear needs to be upgraded to a spatialPixelsDataFrame so that it is a grid # Despite documentation, A SpatialPointsDataFrame object is not a valid input. # first argument(spatialPoints): measurement locations (in our case, the GDD values) # second argument(spatialPixels): where to place the interpolated grid of measurement values. # I need to make an empty grid to hold the data. interpImage <- interpolate(gddThisYear,monarchStaLocs,list(mean=TRUE, variance=TRUE)) message("Latest interp done!") browser() # Assign the result of interpolation as an attribute to the Monarch Station Locations File: 'newcol' cbind(monarchStasWithGDD at data,newcol) } # Write the output file writePointsShape(monarchStasWithGDD,outFileName) } } ############################################################################################################ On 10/26/2010 1:09 PM, Edzer Pebesma wrote:
Are you sure you loaded package intamap before trying this? On 10/26/2010 09:05 PM, Rick Reeves wrote:
Hello:
I note a conflict in the documentation for Intamap/interpolate(),
to know how to resolve:
I have two SpatalPointsDataFrame objects:
growingDegreeObservations<- readShapePoints("OhioGDDStationsGP.shp")
# computed surface phenomena, for individual lat/lon points
OrgStaLocs<- readShapePoints("Org_Locations.shp") #
species observation stations, for DIFFERENT lat/lon points
# Objective: Interpolate the growingDegree points into a SpatialGrid (?)
raster field, and then
# overlay the observation stations on to the field, and assign field
values to each observation station points
# according to intamap documentation, the interpolate() method accepts
two SptatialPointsDataFrame objects, interpolates
# the first ($value field) into a grid, and assigns interpolated
observations to the predicted locations (monarchStaLocs).
interpImage<- interpolate(growingDegreeObservations
,OrgStaLocs,list(mean=TRUE, variance=TRUE))
However, when this line is executed, the error occurs:
Error in function (classes, fdef, mtable) :
unable to find an inherited method for function "interpolate", for
signature "SpatialPointsDataFrame"
I guess that this refers to the first argument, and that the argument
must be a SpatialPixels or SpatialGrid data frame.
Questions:
1) Do I understand the inputs and outputs of interpolate() correctly?
2) How best to transform the incorrect (SpatialPointsDataFrame) argument
to match the required inputs? E.G., create the interpolated
SpatialPixelsDataFrame
of observations before calling interpolate())
3) Is there a better R technique (e.g., direct call to krige()) for
computing my solution?
Thanks, Rick R
Edzer Pebesma Institute for Geoinformatics (ifgi), University of M?nster Weseler Stra?e 253, 48151 M?nster, Germany. Phone: +49 251 8333081, Fax: +49 251 8339763 http://ifgi.uni-muenster.de http://www.52north.org/geostatistics e.pebesma at wwu.de