我的数据结构如下:
winter.dat<-structure(list(ELON = c(-98.02325, -96.66909, -99.33808, -98.70974,
-98.29216, -97.08568, -99.90308, -100.53012, -99.05847, -95.86621,
-97.25452, -102.49713, -96.63121, -97.69394, -96.35404, -94.6244,
-99.64101, -96.81046, -97.26918, -99.27059, -97.0033, -99.34652,
-97.28585, -96.33309, -96.80407, -98.36274, -99.7279, -97.91446,
-95.32596, -95.2487, -95.64617, -94.84896, -95.88553, -96.32027,
-98.03654, -99.80344, -95.65707, -98.49766, -96.71779, -96.42777,
-99.14234, -98.46607, -101.6013, -98.743583, -97.47978, -95.64047,
-96.0024, -98.48151, -99.05283, -96.35595, -99.83331, -101.22547,
-95.54011, -94.8803, -95.45067, -94.78287, -102.8782, -97.76484,
-97.95442, -98.11139, -95.99716, -96.94394, -99.42398, -97.21271,
-99.01109, -95.78096, -97.74577, -98.56936, -94.84437, -97.95553,
-97.60685, -94.82275, -96.91035, -97.20142, -97.95202, -95.60795,
-97.46488, -96.49749, -97.46414, -97.5107, -97.58245, -96.26265,
-95.91473, -97.22924, -96.76986, -97.04831, -95.55976, -95.27138,
-98.96038, -97.15306, -99.36001, -97.58812, -94.79805, -99.0403,
-96.94822, -96.03706, -100.26192, -97.34146, -95.18116, -97.09527,
-96.06982, -96.95048, -94.98671, -95.01152, -99.13755, -96.67895,
-95.22094, -97.52109, -98.52615, -97.98815, -98.77509, -95.1233,
-94.64496, -95.34805, -94.68778, -99.41682, -96.34222), NLAT = c(34.80833,
34.79851, 34.58722, 36.70823, 34.91418, 34.19258, 36.07204, 36.80253,
35.40185, 35.96305, 36.75443, 36.69256, 35.17156, 36.41201, 35.7805,
34.0433, 36.83129, 36.63459, 33.89376, 35.5915, 34.8497, 36.02866,
36.1473, 34.60896, 35.65282, 36.74813, 35.54615, 35.03236, 34.65657,
34.22321, 36.32112, 35.68001, 36.90987, 33.92075, 35.54848, 35.20494,
35.30324, 36.26353, 34.55205, 36.84053, 36.72562, 35.14887, 36.60183,
34.239444, 35.84891, 35.74798, 35.84162, 35.48439, 34.98971,
35.07073, 34.6855, 36.85518, 34.03084, 33.83013, 36.14246, 36.4821,
36.82937, 34.52887, 35.85431, 36.38435, 34.30876, 34.03579, 34.83592,
36.06434, 36.98707, 34.88231, 36.79242, 34.72921, 36.88832, 35.27225,
36.11685, 34.31072, 36.8981, 34.2281, 34.96774, 36.74374, 35.23611,
36.03126, 35.47226, 35.5556, 35.47112, 35.43172, 35.58211, 34.7155,
36.36114, 35.99865, 35.8257, 36.36914, 35.89904, 36.3559, 35.12275,
34.19365, 35.43815, 36.19033, 35.36492, 36.4153, 36.59749, 35.54208,
35.26527, 36.12093, 34.87642, 34.5661, 35.97235, 34.7107, 34.43972,
34.33262, 36.77536, 34.98224, 35.84185, 34.16775, 35.5083, 35.489,
36.011, 34.90092, 34.98426, 36.42329, 36.51806), RAIN_WINTER14 = c(0.7,
1.8, 1.63, 1.14, 1.43, 4.2, 0.76, 1.42, 0.65, 2.42, 0.95, 0.24,
2.82, 1.33, 1.37, 7.5, 1.22, 1.65, 4.3, 0.54, 2.99, 0.78, 1.17,
5.57, 1.85, 0.99, 0.42, 0.69, 5, 4.23, 1.17, 5.82, 1.28, 4.42,
1.22, 0.58, 4.28, 0.85, 2.12, 1.72, 1.41, 0.93, 0.47, 2.28, 1.43,
3.86, 2.69, 1.17, 1.17, 1.6, 1.12, 0.85, 5.27, 7.11, 1.96, 3.12,
0.25, 1.52, 1.19, 1.07, 3.53, 4.95, 0.87, 1.32, 1.26, 4.53, 0.97,
0.47, 2.35, 1.56, 1.22, 7.55, 1.23, 2.98, 0.53, 1.41, 0.61, 1.74,
1.46, 1.62, 1.71, 2.18, 2.43, 1.72, 1.05, 1.39, 2.52, 1.26, 0.61,
1.4, 1.01, 2.13, 4.95, 0.9, 1.34, 1.69, 1.29, 1.56, 4.4, 1.13,
4.82, 2.44, 5.29, 5.68, 1.69, 5.38, 2, 2.54, 1.17, 2.21, 0.67,
4.38, 5.86, 3.79, 6.14, 1.05, 1.03)), .Names = c("ELON", "NLAT",
"RAIN_WINTER14"), row.names = c(NA, -117L), class = "data.frame")
sensor_points<-structure(list(LON = c(-95.91, -97.51, -98.42, -97.51, -97.34,
-97.86, -95.95, -96.09, -96.43, -96.26, -97.11, -98.61, -95.61,
-95.91, -95.91, -94.45, -94.6, -97.51, -95.78, -96.46, -99.5,
-95.78, -98.42, -95.91, -95.94, -94.97, -95.38, -97.86, -97.37,
-97.51, -94.6, -97.51, -97.47, -97.34, -95.78, -97.06, -95.91,
-97.59, -95.66, -97.76, -96.51, -94.66, -99.5, -95.49, -97.22,
-98.24, -95.52, -95.38, -95.26, -96.14), LAT = c(36.12, 35.46,
34.6, 35.46, 35.22, 36.4, 35.63, 35.99, 33.93, 34.13, 34.19,
34.62, 36.31, 36.12, 36.12, 35.33, 35.4, 35.46, 36.03, 36.29,
34.87, 36.03, 34.6, 36.12, 36.73, 35.91, 35.95, 36.4, 35.01,
35.46, 34.89, 35.46, 35.65, 35.22, 36.03, 36.72, 36.12, 35.24,
35.96, 35.51, 34, 35.13, 34.87, 36.28, 34.72, 35.06, 35.47, 35.95,
36.43, 34.01)), .Names = c("LON", "LAT"), row.names = c(NA, 50L
), class = "data.frame")
我正在使用autoKrige()
R 中 automap 包中的函数为我使用以下代码定义的网格生成插值:
ok_state<-map("state","ok",plot=F)
sp1<-seq(min(ok_state$x,na.rm=T),max(ok_state$x,na.rm=T),length=100)
sp2<-seq(min(ok_state$y,na.rm=T),max(ok_state$y,na.rm=T),length=100)
sp<-expand.grid(sp1,sp2)
S<-SpatialPoints(sp)
gridded(S)<-TRUE
proj4string(S)<-"+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
S<-spTransform(S,CRS("+proj=merc +zone=18s +ellps=WGS84 +datum=WGS84"))
要进行插值,我使用以下代码:
coordinates(winter.dat)=~ELON+NLAT
proj4string(winter.dat)="+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
project_winter.dat=spTransform(winter.dat, CRS("+proj=merc +zone=18s +ellps=WGS84 +datum=WGS84"))
kr.grid<-autoKrige(RAIN_WINTER14~1,project_winter.dat,S)
该代码似乎有效。现在,我想使用over()
sp 包中的函数从插值网格的指定单元格中提取预测。我尝试使用以下代码来完成此操作:
coordinates(sensor_points)=~LON+LAT
proj4string(sensor_points)="+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
project_sensor_points=spTransform(sensor_points, CRS("+proj=merc +zone=18s +ellps=WGS84 +datum=WGS84"))
proj4string(kr.grid$krige_output)==proj4string(project_sensor_points)
over(project_sensor_points,kr.grid$krige_output)
但是这段代码产生了一个充满 NA 值的矩阵。这可能是因为kr.grid$krige_output
is aSpatialPointsDataFrame
而不是 a SpatialPixelsDataFrame
or SpatialGridDataFrame
。
有人对如何解决这个问题有任何建议吗?它可以像转换kr.grid$krige_output
为 aSpatialPixelsDataFrame
或 a一样简单SpatialGridDataFrame
。不幸的是,我不知道该怎么做。