如何为R中的预测创建坐标网格 [英] How to create grid for coordinates for Prediction in R
问题描述
我正在使用高斯过程模型进行预测，现在我需要根据数据中的坐标使用Grid文件，但是我没有坐标并且我不知道如何创建它.
我关注了该
#转换为空间像素st_grid<-rasterToPoints(ras，空间= TRUE)gridded(st_grid)<-TRUEst_grid<-as(st_grid，"SpatialPixels") 而且，当我保存网格文件时，如何用网格坐标提示其他列?因为它只显示新的long和lat列
write.csv(st_grid，file ="st_grid.csv") 我不确定您的代码发生了什么，但似乎您的来源设置不正确.我已经更新了上面的代码，以在芝加哥产生网格.我从Google Maps中随机选择了一个起点，并修改了 x_cell 和 y_cell 以在城市范围内生成尺寸合理的地图.
库(sp)图书馆(rgdal)图书馆(光栅)图书馆(传单)图书馆(mapview)站点<-data.frame(lat = c(41.997946，41.960669，41.960669，41.960669,41.909269,41.931841,41.909269,41.910561,41.866129,41.866129)long = c(-87.654561，-87.747456，-87.67459，-87.646438，-87.747456，-87.67459，-87.67459，-87.619112，-87.747456，-87.691617)，站= 1:10)坐标(站)<-〜〜长+纬度#设置投影.它们是经度和纬度，因此请使用WGS84长经投影proj4string(station)<-CRS("+ init = epsg:4326")#使用mapview功能查看车站位置mapview(站)#3.确定原产地#设置原点ori<-SpatialPoints(cbind(-87.872660，41.619136)，proj4string = CRS("+ init = epsg:4326"))#转换ori的投影#使用EPSG:3857(球形墨卡托)ori_t<-spTransform(ori，CRSobj = CRS("+ init = epsg:3857"))#原点四舍五入到最接近的100x_ori<-轮(coordinates(ori_t)[1，1]/100)* 100y_ori<-round(坐标(ori_t)[1，2]/100)* 100#定义x和y轴有多少个像元x_cell<-60y_cell<-80#定义分辨率为1000米单元大小<-1000#创建范围ext<-范围(x_ori，x_ori +(x_cell * cell_size)，y_ori，y_ori +(y_cell * cell_size))#初始化栅格图层ras<-raster(ext)#设置分辨率为res(ras)<-c(cell_size，cell_size)ras []<-0#投影栅格投影(ras)<-CRS("+ init = epsg:3857")#创建交互式地图mapview(station)+ mapview(ras) 这是我最终得到的图像:
关于您的其他问题，我不确定您是否应该将网格与数据结合起来.根据您提到的问题中链接的教程，对于 krige 示例使用数据 meuse 和网格 meuse.grid 作为参数: lzn.kriged<-krige(log(zinc)〜1，meuse，meuse.grid，model = lzn.fit).检查所使用的功能和软件包是否也是这种情况.
如何挑选原产地?此特定代码的原点是网格的左下角，因此我去了Google Maps，并选择了一个稍微超出城市范围的随机点(基于Google的数据)，所以在下方有点偏左不受限制.
I'm using Gaussian Process model for prediction, and I'm now at the point where I need to use Grid file based on the coordinates I have in my data but I don't have one and I don't know how to create it.
I followed the post on this link , but it shows the grid on Pennsylvania not Chicago where my data coordinates located!
So I'm confused which will be the ideal way to create grid file including the other columns in the data.
station <- data.frame(lat = c(41.997946, 41.960669, 41.960669, 41.960669,41.909269,41.931841,41.909269,41.910561,41.866129,41.866129), long = c(-87.654561, -87.747456, -87.67459, -87.646438,-87.747456,-87.67459,-87.67459,-87.619112,-87.747456,-87.691617),station = 1:10)
station
lat long station
1 41.99795 -87.65456 1
2 41.96067 -87.74746 2
3 41.96067 -87.67459 3
4 41.96067 -87.64644 4
5 41.90927 -87.74746 5
6 41.93184 -87.67459 6
7 41.90927 -87.67459 7
8 41.91056 -87.61911 8
9 41.86613 -87.74746 9
10 41.86613 -87.69162 10
The data include more columns such, Hour, Day, Moths, Year, Speed, and these observations are for the 10 locations over 2 months period, but I only put the coordinates here to get an idea how to create the grid.
Here's my steps in creating the grid following the link above:
# Set the projection. They were latitude and longitude, so use WGS84 long-lat projection
proj4string(station) <- CRS("+init=epsg:4326")
# View the station location using the mapview function
mapview(station)
#3. Determine the origin
# Set the origin
ori_t <- SpatialPoints(cbind(-87.67459, 41.99795), proj4string = CRS("+init=epsg:4326"))
# Convert the projection of ori
# Use EPSG: 3857 (Spherical Mercator)
ori_t <- spTransform(ori, CRSobj = CRS("+init=epsg:3857"))
coordinates(ori_t)
#ori_t <- spTransform(ori, CRSobj = CRS("+init=epsg:3857"))
#coordinates(ori_t)
# The origin has been rounded to the nearest 100
x_ori <- round(coordinates(ori_t)[1, 1]/100) * 100
y_ori <- round(coordinates(ori_t)[1, 2]/100) * 100
# Define how many cells for x and y axis
x_cell <- 250
y_cell <- 200
# Define the resolution to be 1000 meters
cell_size <- 1000
# Create the extent
ext <- extent(x_ori, x_ori + (x_cell * cell_size), y_ori, y_ori + (y_cell * cell_size))
# Initialize a raster layer
ras <- raster(ext)
# Set the resolution to be
res(ras) <- c(cell_size, cell_size)
ras[] <- 0
# Project the raster
projection(ras) <- CRS("+init=epsg:3857")
# Create interactive map
mapview(station) + mapview(ras)
But when I view the map, the grid is located in Pennsylvania area not Chicago! Do you have an idea why? I picked for my lat : 41.99795, and my long:-87.67459 , and when I put them on Google map, it shows Chicago area , but not showing the same on R!!
# Convert to spatial pixel
st_grid <- rasterToPoints(ras, spatial = TRUE)
gridded(st_grid) <- TRUE
st_grid <- as(st_grid, "SpatialPixels")
Also, when I save the grid file, how can I clue the other columns with the grid coordinates? Because it's only shows the new long and lat columns
write.csv(st_grid, file = "st_grid.csv")
I'm not sure what happened with your code, but it seems like your origin was set incorrectly. I've updated the code above to produce a grid over Chicago. I've picked a random starting point from Google Maps and modified x_cell and y_cell to produce a reasonably sized map over the city.
library(sp)
library(rgdal)
library(raster)
library(leaflet)
library(mapview)
station <- data.frame(lat = c(41.997946, 41.960669, 41.960669, 41.960669,41.909269,41.931841,41.909269,41.910561,41.866129,41.866129),
long = c(-87.654561, -87.747456, -87.67459, -87.646438,-87.747456,-87.67459,-87.67459,-87.619112,-87.747456,-87.691617),
station = 1:10)
coordinates(station) <- ~long + lat
# Set the projection. They were latitude and longitude, so use WGS84 long-lat projection
proj4string(station) <- CRS("+init=epsg:4326")
# View the station location using the mapview function
mapview(station)
#3. Determine the origin
# Set the origin
ori <- SpatialPoints(cbind(-87.872660, 41.619136), proj4string = CRS("+init=epsg:4326"))
# Convert the projection of ori
# Use EPSG: 3857 (Spherical Mercator)
ori_t <- spTransform(ori, CRSobj = CRS("+init=epsg:3857"))
# The origin has been rounded to the nearest 100
x_ori <- round(coordinates(ori_t)[1, 1]/100) * 100
y_ori <- round(coordinates(ori_t)[1, 2]/100) * 100
# Define how many cells for x and y axis
x_cell <- 60
y_cell <- 80
# Define the resolution to be 1000 meters
cell_size <- 1000
# Create the extent
ext <- extent(x_ori, x_ori + (x_cell * cell_size), y_ori, y_ori + (y_cell * cell_size))
# Initialize a raster layer
ras <- raster(ext)
# Set the resolution to be
res(ras) <- c(cell_size, cell_size)
ras[] <- 0
# Project the raster
projection(ras) <- CRS("+init=epsg:3857")
# Create interactive map
mapview(station) + mapview(ras)
This is the image I get in the end:
As for your other question, I'm not sure you're supposed to combine the grid with your data. According to the tutorial linked in the question you've mentioned, krige for example uses both the data meuse and the grid meuse.grid as arguments: lzn.kriged <- krige(log(zinc) ~ 1, meuse, meuse.grid, model=lzn.fit). Check whether this is the case also for the function and the package you're using.
EDIT: How to pick the origin? The origin in this particular code is the bottom left corner of the grid, so I went to Google Maps and picked a random point that was slightly outside the city limits (based on Google's data), so a bit below and a bit on the left from the limits.
这篇关于如何为R中的预测创建坐标网格的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
