Spatstat, Interpolating data points in R

I have covariate data which is rainfall. I also have a point pattern which represents settlements on a map within a specific area. I have the problem that part of my area of my region contains sea, however there are no settlements on the sea but the rain data does involve the sea area. I wish to set the rain values to those of a nearby point on land.

Any idea how to do this? Maybe to find which rainfall x and y coordinates represent the sea i want to create a data frame that contains the centre points of all the square kilometres and then import the rain values for the locations which match and see which locations have no values. But im stuck. Here is my code:

> window<- data.frame(Lon=c(-1.560367, -1.078330 ), Lat=c( 50.576342, 51.243823))
> coordinates(window) <- ~Lon + Lat
> proj4string(window) = CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0")
> proj4string(winch2) <- latlong
> window <- spTransform(window, bng)
> (floor (coordinates (window) / 1000) + 0.5) * 1000
> W2<- owin(c(431500,464500), c( 75500, 149500))
> Region<-Settlements[W2]     ###Settlements is my data
> rain_im[W2]    ###rain is my covariate as a pixel image
>[W2])  ###Converted this into a pixel image
> `summary(Region)`

Marked planar point pattern:  308 points
Average intensity 1.261261e-07 points per square unit

Coordinates are given to 2 decimal places
i.e. rounded to the nearest multiple of 0.01 units

marks are numeric, of type ‘double’
    Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
   2.375   19.000   47.500  103.029   88.364 5500.000 

Window: rectangle = [431500, 464500] x [75500, 149500] units
Window area = 2.442e+09 square units

1 answer

  • answered 2018-07-18 02:26 Adrian Baddeley

    Please describe your data more clearly. You say you have "covariate data which is rainfall": what format is this data? Is it a pixel image (whose pixel values are rain amounts), or is it a data frame of rain gauge locations and recorded rain amounts? You have some values of the rain covariate at pixels which lie over the sea: what do you want to do with such values? What data do you have that tells you which locations are in the sea and which are on land?

    If you had polygon data defining the coastline, you could convert this to an owin window object, say Land. Then if your rainfall data is a pixel image Rain you could just do ClipRain <- Rain[Land, drop=FALSE] to restrict the pixel image to the land area (assigning NA to any pixel over the sea).

    If you have rainfall data at some locations, say a point pattern RainRecords with marks which are the rainfall values, then Smooth(RainRecords, sigma) is a pixel image containing a kernel smoothed interpolation of rainfall, and nnmark(RainRecords) is a pixel image in which each pixel value is the rainfall at the nearest recorded location. By clipping or sampling these images you can probably get what you need.