2016-03-14 123 views
2

真的希望你能幫助我。預先感謝您對我從這些頁面學到的所有內容。並且道歉,我的技術知識受到兼容性的限制,我試圖學習R中所有的知識。使用大光柵查看錶格(R)

我的目標: 使用(大)柵格查找表反轉。

我有什麼的時刻:

#Observed data (in this case just as a dataframe) 
obs <- data.frame(runif(100,0,1)) 

#Two sets of simulated data (often n >10 000) 
sim.A <- data.frame(runif(1000,0,1)) 

sim.B <- data.frame(runif(1000,0,1)) 

#Calculate the error [cost] for each observed value and every simulated(A) value 
error.fun <- function(x){sqrt((x-sim.A)^2)}    
error <- apply(obs,1,error.fun) 

#Find the position of the min [error] value 
min.func <- function(x){which(x == min(x),arr.ind = F)}  
cost.min <- apply(error,2,min.func) 

#Subset the simulated (B) dataset at the position of the least error[cost.min] 
LUT.values = data.frame(sim.B[cost.min,]) 

我的問題(S):

1)上面的代碼工作從柵格中提取樣本。但是,我需要用整個(ncell> 1Mil)柵格替換採樣的觀測值。我顯然需要優化上述兩個函數(合成一個?),但最接近的我留下了懷疑,因爲結果與採樣數據嘗試相比較差。

我的大柵格嘗試:

#This runs, but I dont think it's working correctly 

crs.UTM <- CRS("+proj=utm +zone=36 +south +datum=WGS84 +units=m +no_defs+ellps=WGS84 +towgs84=0,0,0") 
r <- raster(crs=crs.UTM) 
extent(r) <- extent(0, 100, 0, 100) 
res(r) <- c(1, 1) 
values(r) <- runif(ncell(r), 0, 1) 

#Simulated data (often n >10 000) 
sim.A <- data.frame(runif(1000,0,1)) 

sim.B <- data.frame(runif(1000,0,1)) 

cost.min.func <- function(x){ 
cost <- sqrt((x-sim.A)^2)  
c.min <- sim.B[which(cost == min(cost),arr.ind = FALSE),]} 

LUT.rst <- calc(r,cost.min.func) 

非常感謝

回答

1

我認爲這是你是什麼

library(raster) 
r <- raster(ext=extent(0, 100, 0, 100), res=1, crs="+proj=utm +zone=36 +south +datum=WGS84 +units=m") 

set.seed(0) 
values(r) <- runif(ncell(r), 0, 1) 

sim.A <- runif(1000,0,1) 
sim.B <- runif(1000,0,1) 

cost <- function(x) { 
    y <- abs(x-sim.A) 
    sim.B[which.min(y)] 
} 
x <- calc(r, cost) 

這將需要一段時間在大型數據集上。它應該可以通過使用範圍的x值的第一近似這一點,那麼也許只考慮計算這對於可能具有最小值

+0

謝謝羅伯特幾個細胞。這是完美的,而且要快得多。我的低準確度依然存在,但我會繼續深入研究(不適當的問題等)。儘管感謝這個解決方案。 – Sentinel1b

0

你所缺少的是observed - simulated平均:

rmse <- sqrt(mean((obs-sim)^2))