无法正确使用automap包中的autokrige(R无法很好地读取预测位置)



我正试图使用R来绘制从伊比利亚半岛收集的数据的插值频率图。(像这样https://gis.stackexchange.com/questions/147660/strange-spatial-interpolation-results-from-ordinary-kriging)

我的问题是,由于autokrige函数的atribute new_data中的某种错误,该图没有显示插值数据。

https://cran.r-project.org/web/packages/automap/automap.pdf新数据(_D):包含预测位置的sp对象。new_data可以是一个点集网格或多边形。不得包含NA。如果未提供此对象,则会计算默认值。这是通过取input_data的凸包并在该凸包中放置大约5000个网格单元来完成的

我认为问题是R没有很好地读取转换为多边形的映射,因为如果我避开这个new_data属性,我会得到克里格值的propper图。但我并没有弄清楚伊比利亚半岛的形状。有人能帮帮我吗?我将非常感谢

在这里你可以看到我的数据:http://pastebin.com/QHjn4qjP

实际代码:现在,由于我将数据坐标转换为UTM投影,我没有收到错误消息,但最后一张图没有插值,整个地图显示为单一颜色:(

 setwd("C:/Users/Ruth/Dropbox/R/pruebas")
#Libraries
library(maps)
library(mapdata)
library(automap)
library(sp)
library(maptools)
library(raster)
library(rgdal)
####################MAPA#############
#obtain the grid of the desired country
map<-getData('GADM', country='ESP', level=0)
#convert the grid to projected data
mapa.utm<-spTransform(mapa3,CRSobj =CRS(" +proj=utm +zone=29 +zone=30 +zone=31  +ellps=WGS84"))
###############################Datos#######################
#submit the data
data1<-read.table("FRECUENCIASH.txt",header=T)
head(data1)
attach(data1)
#convert longlat coordinates to UTM
coordinates(data1)<-c("X","Y")
proj4string(data1) = CRS("+proj=longlat +datum=WGS84") 
data1.utm=spTransform(data1, CRS("+proj=utm +zone=29 +zone=30 +zone=31 +ellps=WGS84 "))
######################Kriging interpolation #####################

#Performs an automatic interpolation
kriging_result<-autoKrige(Z~1,data1.utm,mapa.utm,data_variogram = data1.utm) 
#Plot the kriging
result1<-automapPlot(kriging_result$krige_output,"var1.pred",sp.layout = list("sp.points", data1.utm));result1
result2<-plot(kriging_result);result2

您得到的错误与您使用未投影的数据作为自动映射的输入有关。自动地图只能处理投影数据。在谷歌上搜索map projections应该会给你一些背景信息。要将数据投影到合适的投影,可以使用sp包中的spTransform

事实上,它在没有新数据的情况下工作是因为在该对象中没有设置投影,所以自动映射无法警告您。然而,使用latlon输入数据进行自动映射的结果并不可靠。

相关内容

最新更新