[GIS] After I give kring command, this error comes up: data item in gstat object and newdata have different coordinate reference systems

gstatkrigingrsp

> library(gstat)
> library(readxl)
> library(sp)
> library(rgdal)

#2013(FATALITIES)

data <- read.csv("C:/Users/sbs/Desktop/Study/RA/Accident analysis/New folder/Iowa total crashes 2013_cleaned.csv", header = TRUE, sep =",")
plot(data)
attach(data)

> class(data)

[1] "data.frame"

> coordinates(data) <- ~XCOORD+YCOORD

> class(data)

[1] "SpatialPointsDataFrame"
attr(,"package")
[1] "sp"

> variogram1<- variogram(FATALITIES~1, data)
> plot(variogram1)

> fit.variogram1<- fit.variogram(variogram1, vgm(0.002, "Exp"))

> plot(variogram1, fit.variogram1, cutoff = 300000, pch =20, col ="red")

> fit.variogram1

  model       psill    range
1   Exp 0.007896732 13527.15

Until this point everything looks good. My "data" is in CVS format, units meters. Fitted variogram model looks good.

I have made "X" my prediction surface or grid using ArcGIS having 500*500 meter as cell size of fishnet. It is a shapefile.

> x<- readOGR("H:/p/Road_Network_Info/road_fish_lable.shp")

OGR data source with driver: ESRI Shapefile 
Source: "H:\p\Road_Network_Info\road_fish_lable.shp", layer: "road_fish_lable"
with 5006 features
It has 1 fields
Integer64 fields read as strings:  Id 

> plot(x)

> class(x)

[1] "SpatialPointsDataFrame"
attr(,"package")
[1] "sp"

Everything good until now. When I start performing kriging, this error shows up. But the plots of my "data" and Grid (X) overlap each other perfectly.

> krige(FATALITIES~1, data, x, fit.variogram1)

[1] NA
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
Error in predict.gstat(g, newdata = newdata, block = block, nsim = nsim,  : 
  var1 : data item in gstat object and newdata have different coordinate reference systems

I am performing analysis on fatalities in state of Iowa, USA.
Here, both X and data are in meters with same coordinates. i have checked it twice. I am using NAD 1987 Coordinates in ArcGIS to make the required shapefile.

Best Answer

The proj string for your newdata appears to be coordinate reference system 4326, which is in units degrees. If your data is, as you say, in units meters, then they are indeed not in identical coordinate reference systems. Try reprojecting your data to that string returned in the error first and see if that solves the problem:

library(rgdal)

dataReprojected <- spTransform(data, CRS=CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"))

then carry on with your analysis, using dataReprojected.

The error also seems to say the your data has no CRS (it returns two values, one is NA the other a proj string, so maybe you just need to set the CRS of your data in the first place. If your data is in units meters then I don't think that other proj string is correct for your data: find out the CRS of your data and set it...