diff --git a/code/04_buildHeader.Rmd b/code/04_buildHeader.Rmd
index 9e7e932a75b15944bc910182b2c8f073c0238b78..b9874a14e971f7e8bdc91e703c570806f5783c20 100644
--- a/code/04_buildHeader.Rmd
+++ b/code/04_buildHeader.Rmd
@@ -160,7 +160,11 @@ toswap <- c(which(header$Dataset=="Japan" & header$Latitude>90),
header[toswap, c("Latitude", "Longitude")] <- header[toswap, c("Longitude", "Latitude")]
```
-There are `r nrow(header %>% filter(is.na("Location uncertainty (m)")))`. As a first approximation, we assign the median of the respective dataset, as a negative value to indicate this is an estimation, rather than a measure.
+```{r echo=F}
+nouncert <- nrow(header %>% filter(is.na(`Location uncertainty (m)`)))
+```
+
+There are `r nouncert` plots without location uncertainty. As a first approximation, we assign the median of the respective dataset, as a negative value to indicate this is an estimation, rather than a measure.
```{r}
header <- header %>%
left_join(header %>%
@@ -172,7 +176,11 @@ header <- header %>%
`Location uncertainty (m)`)) %>%
dplyr::select(-loc.uncer.median)
```
-There are still `r nrow(header %>% filter(is.na("Location uncertainty (m)")))` plots with no estimation of location uncertainty.
+
+```{r echo=F}
+nouncert <- nrow(header %>% filter(is.na(`Location uncertainty (m)`)))
+```
+There are still `r nouncert` plots with no estimation of location uncertainty.
\newline
Assign plot size to plots in the Patagonia dataset (input of Ana Cingolani)
```{r}
@@ -254,8 +262,8 @@ header <- header %>%
)) %>%
mutate(`Plants recorded`=factor(`Plants recorded`, exclude = "#N/A"))
```
-Align labels consortium
-```{r message=F}
+Align consortium labels to those in sPlot's consortium archive
+```{r message=F, warning=F}
databases <- read_csv("/data/sPlot/users/Francesco/_sPlot_Management/Consortium/Databases.out.csv")
header <- header %>%
@@ -358,6 +366,15 @@ header <- header %>%
labels=c("AF", "AN", "AU", "EU", "NA", "SA"))) %>%
dplyr::select(-continent)
```
+Summarize
+```{r message=F, echo=F}
+knitr::kable(header %>%
+ group_by(CONTINENT) %>%
+ summarize(num.plot=n()),
+ caption="Number of plots per continent") %>%
+ kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
+ full_width = F, position = "center")
+```
## 4.2 Assign to sBiomes
Performed in EVE HPC cluster using function `A98_PredictorsExtract.R`. Divided in 99 chunks.
@@ -398,8 +415,17 @@ header <- header %>%
```
-There are `r sum(is.na(sBiomes.out$Name))` unassigned plots.
-
+There are `r sum(is.na(sBiomes.out$Name))` unassigned plots.
+
+Summarize:
+```{r message=F, echo=F}
+knitr::kable(header %>%
+ group_by(sBiome) %>%
+ summarize(num.plot=n()),
+ caption="Number of plots per Biome") %>%
+ kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
+ full_width = F, position = "center")
+```
## 4.3 Extract elevation
@@ -418,6 +444,7 @@ header.tiles <- header %>%
```
There are `r nrow(header.tiles)` plots out of `r nrow(header)` plots with Location uncertainty <= 50km (or absent). The total number of tiles is `r nlevels(header.tiles$tilenam)`.
+Performed in EVE HPC cluster using function `A97_ElevationExtract.R`. Divided in 99 chunks.
```{r eval=F}
cl <- makeForkCluster(5, outfile="")
registerDoParallel(cl)
@@ -427,8 +454,7 @@ clusterEvalQ(cl, {
library(raster)
library(sp)
library(elevatr)
- library(dplyr)
- })
+ library(dplyr)})
# Divided in 99 chunks
elevation.out <- foreach(i=1:99, .combine=rbind) %dopar% {
@@ -465,11 +491,9 @@ clusterEvalQ(cl, {
library(raster)
library(sp)
library(elevatr)
- library(dplyr)
- })
+ library(dplyr)})
#Loop over all plots
-#for(i in nrow(sp.tile0)){
elevation.failed <- foreach(i=1:nrow(sp.tile0), .packages=c('raster'), .combine=rbind) %dopar% {
sp.tile <- sp.tile0[i,]
tryCatch(raster.tile <- get_elev_raster(sp.tile, z=10,
@@ -493,20 +517,16 @@ elevation.failed <- foreach(i=1:nrow(sp.tile0), .packages=c('raster'), .combine=
elev.tile <- raster::extract(raster.tile, sp.tile, small=T)
elev.tile.buffer <- raster::extract(raster.tile, sp.tile,
buffer=sp.tile$`Location uncertainty (m)`, small=T)
- #elev.q95 <- t(round(sapply(elev.tile.buffer,stats::quantile, probs=c(0.025, 0.5, 0.975), na.rm=T)))
elev.q95 <- t(round(mapply( quantile,
x=elev.tile.buffer,
- #center=elev.tile,
- probs=rep(c(0.025, 0.5, 0.975), each=length(elev.tile)),
- #loc.uncert=sp.tile$`Location uncertainty (m)`,
- na.rm=T)))
- output.tile <- data.frame(PlotObservationID=sp.tile$PlotObservationID,
+ probs=rep(c(0.025, 0.5, 0.975), each=length(elev.tile)), na.rm=T)))
+ output.tile <- data.frame(PlotObservationID=sp.tile$PlotObservationID,
elevation=round(elev.tile),
elev.q95,
DEM.res=res(raster.tile)[1]) %>%
- rename(Elevation_q2.5=X2.5., Elevation_median=X50., Elevation_q97.5=X97.5.)
- return(output.tile)
- }
+ rename(Elevation_q2.5=X2.5., Elevation_median=X50., Elevation_q97.5=X97.5.)
+ return(output.tile)
+ }
}
stopCluster(cl)
save(elevation.failed, file = "../_derived/elevatr/elevation_missing.RData")
@@ -548,25 +568,32 @@ write_csv(elevation.out, path ="../_derived/elevatr/elevation.out.csv")
Reimport output, attach to header and check
-```{r message=F}
+```{r message=F, echo=F}
elevation.out <- read_csv("../_derived/elevatr/elevation.out.csv")
knitr::kable(head(elevation.out,10),
caption="Example of elevation output (only first 10 shown)") %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = F, position = "center")
-summary(elevation.out)
```
-There are `r sum(is.na(elevation.out$elevation))` plots without elevation info, corresponding to `r round(sum(is.na(elevation.out$elevation))/nrow(header)*100,1)`% of total.
-There are `r sum(elevation.out$elevation < -1, na.rm=T)` plots with elevation below sea level.
+There are `r sum(is.na(elevation.out$Elevation_median))` plots without elevation info, corresponding to `r round(sum(is.na(elevation.out$Elevation_median))/nrow(header)*100,1)`% of total.
+There are `r sum(elevation.out$Elevation_median < -1, na.rm=T)` plots with elevation below sea level.
+\newline
+Join elevation data (only median)
+```{r}
+header <- header %>%
+ left_join(elevation.out %>%
+ dplyr::select(PlotObservationID, Elevation_median) %>%
+ rename(elevation_dem=Elevation_median),
+ by="PlotObservationID")
-```{r message=F}
+```
+Summary and check
+```{r message=F, echo=F}
elevation.out <- read_csv("../_derived/elevatr/elevation.out.csv")
knitr::kable(header %>%
- dplyr::select(PlotObservationID, Dataset, `GIVD ID`) %>%
- left_join(elevation.out,
- by="PlotObservationID") %>%
- filter(elevation< -1) %>%
+ dplyr::select(PlotObservationID, elevation_dem, Dataset, `GIVD ID`) %>%
+ filter(elevation_dem < -1) %>%
group_by( `GIVD ID`, Dataset) %>%
summarize(num.plot=n()) %>%
ungroup() %>%
@@ -580,7 +607,7 @@ knitr::kable(header %>%
```{r}
summary(elevation.out %>%
- dplyr::select(-PlotObservationID))
+ dplyr::select(-PlotObservationID, -elevation))
```
@@ -603,6 +630,7 @@ ggplot(data=mydata) +
```
+
# 5 Map of plots
```{r, message=F, warning=F}
@@ -635,7 +663,7 @@ w3a <- ggplot() +
```
Graph of plot density (hexagons)
-```{r, fig.align="center", fig.width=8, fig.height=4, cache=T}
+```{r, message=F, fig.align="center", fig.width=8, fig.height=4, cache=T}
header2 <- header %>%
filter(!is.na(Longitude) | !is.na(Latitude)) %>%
dplyr::select(PlotObservationID, Latitude, Longitude, `GIVD ID`) %>%
@@ -672,7 +700,7 @@ ggsave(filename="../_pics/PlotDensityLog10_vir.png", width = 15, height = 7, uni
```
Graph of plot location by Dataset
-```{r, fig.align="center", fig.width=8, fig.height=4, cache=T}
+```{r, fig.align="center", fig.width=8, fig.height=4, cache=T, message=F}
header.sf <- header.shp %>%
st_as_sf() %>%
st_transform(crs = "+proj=eck4")
@@ -708,6 +736,12 @@ save(header.out, file = "../_output/header_splot3.0")
```
+```{r message=F, echo=F}
+knitr::kable(header.out %>% slice(1:20),
+ caption="Example of header (first 20 rows)") %>%
+ kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
+ full_width = F, position = "center")
+```
## Supplementary Material
### ANNEX 1 - Ancillary function - `PredExtr`
@@ -895,7 +929,10 @@ ElevationExtract <- function(header, output, ncores, chunk.i){
stopCluster(cl)
}
-
+```
+### ANNEX 3 - SessionInfo()
+```{r}
+sessionInfo()
```