Skip to content
Snippets Groups Projects
Commit ae8bd4fd authored by Francesco Sabatini's avatar Francesco Sabatini
Browse files

refined buildHeader

parent d707dc27
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 70 additions and 33 deletions
code/04_buildHeader.Rmd
+ 70
− 33
View file @ ae8bd4fd
......@@ -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()
```
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment