set.seed(123)
require(sp); require(sf)
data(meuse)
data(meuse.grid)
### Data
y <- log(meuse[,"zinc"])
coords <- meuse[,c("x","y")]
x <- data.frame(dist = meuse[,"dist"],
ffreq2 = as.integer(meuse$ffreq == 2),
ffreq3 = as.integer(meuse$ffreq == 3))
### Data at prediction sites
coords0 <- meuse.grid[,c("x","y")]
x0 <- data.frame(dist = meuse.grid[,"dist"],
ffreq2 = as.integer(meuse.grid$ffreq == 2),
ffreq3 = as.integer(meuse.grid$ffreq == 3))
### Holdout validation optimizing the number of spatial scales
mod_hv <- cf_lm_hv(y = y, x = x, coords = coords, add_learn = "none")
### Spatial modeling and prediction
mod <- cf_lm(y = y, x = x, x0 = x0, coords = coords, coords0 = coords0,
mod_hv = mod_hv)
mod
### Mapping predictive mean and standard deviations (SD)
meuse.grid$pred <- mod$pred0$pred
meuse.grid$pred_sd<- mod$pred0$pred_sd
meuse.grid_sf <- st_as_sf(meuse.grid, coords = c("x","y"))
plot(meuse.grid_sf[,"pred"], pch = 15, cex = 0.5, nbreaks = 20) # Predictive mean
plot(meuse.grid_sf[,"pred_sd"], pch = 15, cex = 0.5, nbreaks = 20)# Predictive SD
### Multiscale spatial pattern/feature extraction
mod_s1<- sp_scalewise(mod,bw_range=c(1000,Inf)) # Large scale (1000 <= bandwdith)
mod_s2<- sp_scalewise(mod,bw_range=c(500,1000)) # Middle scale (500 <= bandwdith <= 1000)
mod_s3<- sp_scalewise(mod,bw_range=c(0,500)) # Small scale (bandwdith <= 500)
z1 <- mod_s1$pred0$pred # Predictive mean
z2 <- mod_s2$pred0$pred
z3 <- mod_s3$pred0$pred
z1_sd <- mod_s1$pred0$pred_sd # Predictive SD
z2_sd <- mod_s2$pred0$pred_sd
z3_sd <- mod_s3$pred0$pred_sd
meuse.grid_sf3 <- cbind(meuse.grid_sf, z1, z2, z3, z1_sd, z2_sd, z3_sd)
plot(meuse.grid_sf3[,c("z1","z2","z3")], pch = 15,
cex = 0.5, nbreaks = 20,key.pos=4,axes=TRUE) # Predictive means
plot(meuse.grid_sf3[,c("z1_sd","z2_sd","z3_sd")], pch = 15,
cex = 0.5, nbreaks = 20,key.pos=4,axes=TRUE) # Predictive SD
Run the code above in your browser using DataLab