if (FALSE) {
library(sf)
library(terra)
library(caret)
library(viridis)
# prepare sample data:
data(cookfarm)
dat <- aggregate(cookfarm[,c("VW","Easting","Northing")],
by=list(as.character(cookfarm$SOURCEID)),mean)
pts <- st_as_sf(dat,coords=c("Easting","Northing"),crs=26911)
pts$ID <- 1:nrow(pts)
set.seed(100)
pts <- pts[1:30,]
studyArea <- rast(system.file("extdata","predictors_2012-03-25.tif",package="CAST"))[[1:8]]
trainDat <- extract(studyArea,pts,na.rm=FALSE)
trainDat <- merge(trainDat,pts,by.x="ID",by.y="ID")
# visualize data spatially:
plot(studyArea)
plot(studyArea$DEM)
plot(pts[,1],add=TRUE,col="black")
# train a model:
set.seed(100)
variables <- c("DEM","NDRE.Sd","TWI")
model <- train(trainDat[,which(names(trainDat)%in%variables)],
trainDat$VW, method="rf", importance=TRUE, tuneLength=1,
trControl=trainControl(method="cv",number=5,savePredictions=T))
print(model) #note that this is a quite poor prediction model
prediction <- predict(studyArea,model,na.rm=TRUE)
plot(varImp(model,scale=FALSE))
#...then calculate the AOA of the trained model for the study area:
AOA <- aoa(studyArea, model)
plot(AOA)
plot(AOA$AOA)
#... or if preferred calculate the aoa and the LPD of the study area:
AOA <- aoa(studyArea, model, LPD = TRUE, maxLPD = 1)
plot(AOA$LPD)
####
#The AOA can also be calculated without a trained model.
#All variables are weighted equally in this case:
####
AOA <- aoa(studyArea,train=trainDat,variables=variables)
####
# The AOA can also be used for models trained via mlr3 (parameters have to be assigned manually):
####
library(mlr3)
library(mlr3learners)
library(mlr3spatial)
library(mlr3spatiotempcv)
library(mlr3extralearners)
# initiate and train model:
train_df <- trainDat[, c("DEM","NDRE.Sd","TWI", "VW")]
backend <- as_data_backend(train_df)
task <- as_task_regr(backend, target = "VW")
lrn <- lrn("regr.randomForest", importance = "mse")
lrn$train(task)
# cross-validation folds
rsmp_cv <- rsmp("cv", folds = 5L)$instantiate(task)
## predict:
prediction <- predict(studyArea,lrn$model,na.rm=TRUE)
### Estimate AOA
AOA <- aoa(studyArea,
train = as.data.frame(task$data()),
variables = task$feature_names,
weight = data.frame(t(lrn$importance())),
CVtest = rsmp_cv$instance[order(row_id)]$fold)
}
Run the code above in your browser using DataLab