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(gddta...@data$yr)
numYears <- length(theYears)
for (curYear in theYears)
{
gddThisYear <- gddTable[which(gddta...@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(monarchstaswith...@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
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo