Skip to content

Commit

Permalink
harmonize with main
Browse files Browse the repository at this point in the history
Merge branch 'main' into spcv_fixes

# Conflicts:
#	geocompr.bib
  • Loading branch information
jannes-m committed Apr 20, 2022
2 parents c1f763b + aa4039a commit 631b7c4
Show file tree
Hide file tree
Showing 3 changed files with 24 additions and 20 deletions.
30 changes: 15 additions & 15 deletions 15-eco.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,7 @@ To do so, we will bring together concepts presented in previous chapters and eve
Fog oases are one of the most fascinating vegetation formations we have ever encountered.
These formations, locally termed *lomas*, develop on mountains along the coastal deserts of Peru and Chile.^[Similar vegetation formations develop also in other parts of the world, e.g., in Namibia and along the coasts of Yemen and Oman [@galletti_land_2016].]
The deserts' extreme conditions and remoteness provide the habitat for a unique ecosystem, including species endemic to the fog oases.
Despite the arid conditions and low levels of precipitation of around 30-50 mm per year on average, fog deposition increases the amount of water available to plants during austal winter.
Despite the arid conditions and low levels of precipitation of around 30-50 mm per year on average, fog deposition increases the amount of water available to plants during austral winter.
This results in green southern-facing mountain slopes along the coastal strip of Peru (Figure \@ref(fig:study-area-mongon)).
This fog, which develops below the temperature inversion caused by the cold Humboldt current in austral winter, provides the name for this habitat.
Every few years, the El Niño phenomenon brings torrential rainfall to this sun-baked environment [@dillon_lomas_2003].
Expand Down Expand Up @@ -424,7 +424,7 @@ text(tree_mo, pretty = 0)
dev.off()
```

```{r tree, echo=FALSE, fig.cap="Simple example of a decision tree with three internal nodes and four terminal nodes.", fig.scap="Simple example of a decision tree."}
```{r tree, echo=FALSE, fig.cap="Simple example of a decision tree with three internal nodes and four terminal nodes.", out.width="60%", fig.scap="Simple example of a decision tree."}
knitr::include_graphics("figures/15_tree.png")
```

Expand Down Expand Up @@ -475,7 +475,7 @@ Therefore, we tune the hyperparameters\index{hyperparameter} for a good spatial
<!-- If we used more than one repetition (say 2) we would retrieve multiple optimal tuned hyperparameter combinations (say 2) -->

Having already constructed the input variables (`rp`), we are all set for specifying the **mlr3**\index{mlr3 (package)} building blocks (task, learner, and resampling).
For specifying a spatial task, we use again the **mlr3spatiotempcv** [@schratz_mlr3spatiotempcv_2021] package (section \@ref(spatial-cv-with-mlr3)).
For specifying a spatial task, we use again the **mlr3spatiotempcv** package [@schratz_mlr3spatiotempcv_2021 & Section \@ref(spatial-cv-with-mlr3)].
Since our response (`sc`) is numeric, we use a regression\index{regression} task.

```{r 15-eco-20, eval=FALSE}
Expand All @@ -486,7 +486,7 @@ task = mlr3spatiotempcv::TaskRegrST$new(

Using an `sf` object as the backend automatically provides the geometry information needed for the spatial partitioning later on.
Additionally, we got rid of the columns `id` and `spri` since these variables should not be used as predictors in the modeling.
Next, we go on to contruct the a random forest\index{random forest} learner from the **ranger** package.
Next, we go on to construct the a random forest\index{random forest} learner from the **ranger** package.

```{r 15-eco-21, eval=FALSE}
lrn_rf = lrn("regr.ranger", predict_type = "response")
Expand Down Expand Up @@ -519,11 +519,11 @@ search_space = paradox::ps(
Having defined the search space, we are all set for specifying our tuning via the `AutoTuner()` function.
Since we deal with geographic data, we will again make use of spatial cross-validation to tune the hyperparameters\index{hyperparameter} (see Sections \@ref(intro-cv) and \@ref(spatial-cv-with-mlr3)).
Specifically, we will use a five-fold spatial partitioning with only one repetition (`rsmp()`).
In each of these spatial partitions, we run 50 models (`trm()`) while using randomly selected hyperparameter configurations (`tnr`) within predefined limits (`seach_space`) to find the optimal hyperparameter\index{hyperparameter} combination.
In each of these spatial partitions, we run 50 models (`trm()`) while using randomly selected hyperparameter configurations (`tnr()`) within predefined limits (`seach_space`) to find the optimal hyperparameter\index{hyperparameter} combination [see also Section \@ref(svm) and https://mlr3book.mlr-org.com/optimization.html#autotuner, @becker_mlr3_2022].
The performance measure is the root mean squared error (RMSE\index{RMSE}).

```{r 15-eco-22, eval=FALSE}
at = mlr3tuning::AutoTuner$new(
autotuner_rf = mlr3tuning::AutoTuner$new(
learner = lrn_rf,
# spatial partitioning
resampling = mlr3::rsmp("spcv_coords", folds = 5),
Expand All @@ -541,9 +541,9 @@ at = mlr3tuning::AutoTuner$new(
Calling the `train()`-method of the `AutoTuner`-object finally runs the hyperparameter\index{hyperparameter} tuning, and will find the optimal hyperparameter\index{hyperparameter} combination for the specified parameters.

```{r 14-eco-24, eval=FALSE}
# hyperparamter tuning
# hyperparameter tuning
set.seed(0412022)
at$train(task)
autotuner_rf$train(task)
#>...
#> INFO [11:39:31.375] [mlr3] Finished benchmark
#> INFO [11:39:31.427] [bbotk] Result of batch 50:
Expand All @@ -558,16 +558,16 @@ at$train(task)
```

```{r 15-eco-25, eval=FALSE, echo=FALSE}
saveRDS(at, "extdata/15-tune.rds")
saveRDS(autotuner_rf, "extdata/15-tune.rds")
```

```{r 15-eco-26, echo=FALSE, eval=FALSE}
tune = readRDS("extdata/15-tune.rds")
autotuner_rf = readRDS("extdata/15-tune.rds")
```

An `mtry` of 4, a `sample.fraction` of 0.9, and a `min.node.size` of 7 represent the best hyperparameter\index{hyperparameter} combination.
An RMSE\index{RMSE} of
<!-- `r # round(tune$tuning_result$regr.rmse, 2)` -->
<!-- `r # round(autotuner_rf$tuning_result$regr.rmse, 2)` -->
0.38
is relatively good when considering the range of the response variable which is
<!-- `r # round(diff(range(rp$sc)), 2)` -->
Expand All @@ -581,17 +581,17 @@ To do so, we only need to run the `predict` method of our fitted `AutoTuner` obj

```{r 15-eco-27, eval=FALSE}
# predicting using the best hyperparameter combination
at$predict(task)
autotuner_rf$predict(task)
```

The `predict` method will apply the model to all observations used in the modeling.
Given a multilayer `SpatRaster` containing rasters named as the predictors used in the modeling, `terra::predict()` will also make spatial predictions, i.e., predict to new data.

```{r 15-eco-28, eval=FALSE}
pred = terra::predict(ep, model = at, fun = predict)
pred = terra::predict(ep, model = autotuner_rf, fun = predict)
```

```{r rf-pred, echo=FALSE, fig.cap="Predictive mapping of the floristic gradient clearly revealing distinct vegetation belts.", fig.width = 10, fig.height = 10, fig.scap="Predictive mapping of the floristic gradient."}
```{r rf-pred, echo=FALSE, fig.cap="Predictive mapping of the floristic gradient clearly revealing distinct vegetation belts.", out.width="60%", fig.scap="Predictive mapping of the floristic gradient."}
# # restrict the prediction to your study area
# pred = terra::mask(pred, terra::vect(study_area)) |>
# terra::trim()
Expand Down Expand Up @@ -632,7 +632,7 @@ newdata = as.data.frame(as.matrix(ep))
colSums(is.na(newdata)) # 0 NAs
# but assuming there were 0s results in a more generic approach
ind = rowSums(is.na(newdata)) == 0
tmp = at$predict_newdata(newdata = newdata[ind, ], task = task)
tmp = autotuner_rf$predict_newdata(newdata = newdata[ind, ], task = task)
newdata[ind, "pred"] = data.table::as.data.table(tmp)[["response"]]
pred_2 = ep$dem
# now fill the raster with the predicted values
Expand Down
12 changes: 8 additions & 4 deletions _15-ex.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@ E1. Run a NMDS\index{NMDS} using the percentage data of the community matrix.
Report the stress value and compare it to the stress value as retrieved from the NMDS using presence-absence data.
What might explain the observed difference?

```{r 15-ex-e1, eval=FALSE}
```{r 15-ex-e1, message=FALSE}
data("comm", package = "spDataLarge")
pa = decostand(comm, "pa")
pa = pa[rowSums(pa) != 0, ]
Expand All @@ -35,6 +35,7 @@ nmds_pa$stress
nmds_per$stress
```

```{asis, message=FALSE}
The NMDS using the presence-absence values yields a better result (`nmds_pa$stress`) than the one using percentage data (`nmds_per$stress`).
This might seem surprising at first sight.
On the other hand, the percentage matrix contains both more information and more noise.
Expand All @@ -47,6 +48,7 @@ The point here is that percentage data as specified during a field campaign migh
This again introduces noise which in turn will worsen the ordination result.
Still, it is a valuable information if one species had a higher frequency or coverage in one plot than another compared to just presence-absence data.
One compromise would be to use a categorical scale such as the Londo scale.
```

E2. Compute all the predictor rasters\index{raster} we have used in the chapter (catchment slope, catchment area), and put them into a `SpatRaster`-object.
Add `dem` and `ndvi` to it.
Expand All @@ -55,7 +57,7 @@ Finally, construct a response-predictor matrix.
The scores of the first NMDS\index{NMDS} axis (which were the result when using the presence-absence community matrix) rotated in accordance with elevation represent the response variable, and should be joined to `random_points` (use an inner join).
To complete the response-predictor matrix, extract the values of the environmental predictor raster object to `random_points`.

```{r 15-ex-e2, eval=FALSE}
```{r 15-ex-e2}
# first compute the terrain attributes we have also used in the chapter
library(dplyr)
library(terra)
Expand Down Expand Up @@ -118,7 +120,7 @@ Parallelize\index{parallelization} the tuning level (see Section \@ref(svm)).
Report the mean RMSE\index{RMSE} and use a boxplot to visualize all retrieved RMSEs.
Please not that this exercise is best solved using the mlr3 functions `benchmark_grid()` and `benchmark()` (see https://mlr3book.mlr-org.com/perf-eval-cmp.html#benchmarking for more information).

```{r 15-ex-e3, eval=FALSE}
```{r 15-ex-e3, message=FALSE}
library(dplyr)
library(future)
library(mlr3)
Expand Down Expand Up @@ -190,7 +192,7 @@ tictoc::toc()
# stop parallelization
future:::ClusterRegistry("stop")
# save your result, e.g. to
saveRDS(bmr, file = "extdata/15_bmr.rds")
# saveRDS(bmr, file = "extdata/15_bmr.rds")
# mean RMSE
bmr$aggregate(measures = msr("regr.rmse"))
Expand All @@ -213,5 +215,7 @@ ggplot(data = d, mapping = aes(x = learner_id, y = regr.rmse)) +
labs(y = "RMSE", x = "model")
```

```{asis, message=FALSE}
In fact, `lm` performs at least as good the random forest model, and thus should be preferred since it is much easier to understand and computationally much less demanding (no need for fitting hyperparameters).
But keep in mind that the used dataset is small in terms of observations and predictors and that the response-predictor relationships are also relatively linear.
```
2 changes: 1 addition & 1 deletion geocompr.bib
Original file line number Diff line number Diff line change
Expand Up @@ -2193,4 +2193,4 @@ @article{lang_mlr3_2019
pages = {1903},
journal = {Journal of Open Source Software},
number = {44}
}
}

0 comments on commit 631b7c4

Please sign in to comment.