In R, and RStudio, I have run his code:
#ordinary kriging
mydata.krig <- krige(Clay~1, train_data, newdata=mask, vgm)
names(mydata.krig)
names(mydata.krig)[1] <- "Clay.pred"
names(mydata.krig)
min(mydata.krig$Clay.pred, na.rm=T); max(mydata.krig$Clay.pred, na.rm=T)
ggplot() +
geom_stars(data = mydata.krig["Clay.pred"]) +
scale_fill_gradient(low = "yellow", high = "dark blue", limits = c(15.5486,61.41131)) +
geom_sf(data = thessaly_smus, color = "black", fill = NA, size = 1) +
labs(title = "%Clay Predicted Values - Ordinary Kriging, OK")
With ggplot
, I have also acquired a map of the predicted values of %Clay
. I want to extract this map as a raster file from R in order to open it in a GIS, but only the area of interest (coloured polygons). Can you provide me accurate code for this?
Updated* With the help of @I_O, I wrote the code as you can see in comments, but again this is not my outcome that I want to take.I took the following outcome:
But I want to create a raster file like this that will contain the clay predictions:
mydata.krig
#> stars object with 2 dimensions and 2 attributes
#> attribute(s):
#> Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
#> Clay.pred 15.80128 30.74017 36.38344 36.62676 43.28872 61.4972 63587
#> var1.var 85.52382 92.99561 96.14868 101.20852 100.43495 193.6164 63587
#> dimension(s):
#> from to offset delta refsys x/y
#> x 1 346 267586 487.021 GGRS87 / Greek Grid [x]
#> y 1 248 4442456 -487.271 GGRS87 / Greek Grid [y]
as mydata.krieg is a stars
object, it’s quite easy to convert it to raster
raster_ext <- sf::st_bbox(mydata.krig)
r <- raster::raster(t(mydata.krig[[1]]),
xmn = raster_ext[1], xmx = raster_ext[3],
ymn = raster_ext[2], ymx=raster_ext[4],
crs = sf::st_crs(mydata.krig)$proj4string)
terra::plot(r, col = grDevices::gray.colors(40))
and write it to file
raster::writeRaster(r, filename = "clay.tiff", overwrite = TRUE)