Skip to content

Commit

Permalink
add pkgdown, add ensemble logic #70
Browse files Browse the repository at this point in the history
  • Loading branch information
mikejohnson51 committed Jul 6, 2023
1 parent 204a23c commit 9aedaf3
Show file tree
Hide file tree
Showing 14 changed files with 392 additions and 578 deletions.
4 changes: 4 additions & 0 deletions .Rbuildignore
Original file line number Diff line number Diff line change
Expand Up @@ -12,3 +12,7 @@ help
^\.github$
.lazytest
img
^_pkgdown\.yml$
^docs$
^pkgdown$
^vignettes/articles$
1 change: 1 addition & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -6,3 +6,4 @@ sampleFiles
inProgress
data_backup.DS_Store
.lazytest
docs
8 changes: 6 additions & 2 deletions DESCRIPTION
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,10 @@ Type: Package
Title: climateR
Description: A set of functions to subset climate data.
Version: 0.3.1
Authors@R: person("Mike", "Johnson", role = c("aut", "cre"), email = "mikecp11@gmail.com")
Authors@R: c(person("Mike", "Johnson",
role = c("aut", "cre"),
email = "mikecp11@gmail.com"),
person("ESIP", role = "fnd"))
Maintainer: Mike Johnson <mikecp11@gmail.com>
BugReports: https://github.com/mikejohnson51/climateR/issues
URL: https://github.com/mikejohnson51/climateR
Expand All @@ -27,7 +30,8 @@ Suggests:
httr,
AOI,
covr,
testthat (>= 3.0.0)
testthat (>= 3.0.0),
rmarkdown
Remotes:
mikejohnson51/AOI
Config/testthat/edition: 3
12 changes: 8 additions & 4 deletions R/climater_filter.R
Original file line number Diff line number Diff line change
Expand Up @@ -123,10 +123,6 @@ climater_filter <- function(id = NULL,
scenario <- c("historical", scenario)
}

# if ("total" %in% catalog$scenario) {
# scenario <- c("total", scenario)
# }

if(!is.null(scenario)){
catalog <- filter(catalog, scenario %in% !!scenario)
}
Expand All @@ -135,6 +131,13 @@ climater_filter <- function(id = NULL,

### 3 ---- ensemble filter

if(length(ensemble) != length(model)){
catalog = catalog %>%
group_by(model, ensemble) %>%
slice(1) %>%
ungroup()
} else {

u <- unique(catalog$ensemble)

if (length(u) > 1 & is.null(ensemble)) {
Expand All @@ -156,6 +159,7 @@ climater_filter <- function(id = NULL,
)
}
}
}



Expand Down
279 changes: 11 additions & 268 deletions README.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,6 @@ library(climateR)
[![Dependencies](https://img.shields.io/badge/dependencies-7/25-orange?style=flat)](#)
[![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](https://choosealicense.com/licenses/mit/)
[![Project Status: Active](https://www.repostatus.org/badges/latest/active.svg)](https://www.repostatus.org/#active)
[![LifeCycle](https://img.shields.io/badge/lifecycle-experimental-orange)](https://lifecycle.r-lib.org/articles/stages.html#experimental)
[![codecov](https://codecov.io/gh/mikejohnson51/climateR/branch/master/graph/badge.svg?token=7zs6C91SDw)](https://codecov.io/gh/mikejohnson51/climateR)
<!-- badges: end -->

Expand All @@ -55,283 +54,27 @@ remotes::install_github("mikejohnson51/AOI") # suggested!
remotes::install_github("mikejohnson51/climateR")
```

# Usful Packages for climate data
# Basic Usage

## Rainfall from GridMet for Colorado
### October 29,1991 - November 6, 1991

```{r}
library(AOI)
library(climateR)
library(tidyterra)
library(ggplot2)
library(terra)
library(tidyr)
library(sf)
```

# Examples

The climateR package is supplemented by the [AOI](https://github.com/mikejohnson51/AOI) framework established in the `AOI` R package.

To get a climate product, an area of interest must be defined:
library(climateR)
```{r}
AOI = aoi_get(state = "NC")
AOI = aoi_get(state = "CO", county = "all")
plot(AOI$geometry)
```

Here we are loading a polygon for the state of North Carolina More examples of constructing AOI calls can be found [here](https://mikejohnson51.github.io/AOI/).

With an AOI, we can construct a call to a dataset for a parameter(s) and date(s) of choice. Here we are querying the PRISM dataset for maximum and minimum temperature on October 29, 2018:

```{r}
system.time({
p = getPRISM(AOI, varname = c('tmax','tmin'), startDate = "2018-10-29")
})
```

```{r, echo = FALSE}
ggplot() +
geom_spatraster(data = rast(p)) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white"
) + theme_minimal()
```

# Data from known bounding coordinates

`climateR` offers support for `sf`, `sfc`, and `bbox` objects. Here we are requesting wind velocity data for the four corners region of the USA by bounding coordinates.

```{r}
AOI = st_as_sfc(st_bbox(c(xmin = -112, xmax = -105, ymax = 39, ymin = 34), crs = 4326))
g = getGridMET(AOI,
varname = "vs",
startDate = "2018-09-01")
ggplot() +
geom_spatraster(data = rast(g)) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white"
) +
geom_sf(data = AOI::aoi_get(state = c("CO", "NM", "AZ", "UT")), fill = NA) +
theme_minimal()
```

# Data through time ...

In addition to multiple variables we can request variables through time, here let's look at the gridMET rainfall for the Gulf Coast during Hurricane Harvey:

```{r, fig.width= 15}
harvey = getGridMET(aoi_get(state = c("TX", "FL")),
varname = "pr",
startDate = "2017-08-20", endDate = "2017-08-31")
ggplot() +
geom_spatraster(data = harvey$precipitation_amount) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white") +
theme_minimal()
```

# Climate Projections

Some sources are downscaled Global Climate Models (GCMs). These allow you to query forecasted ensemble members from different models and/or climate scenarios. One example is from the MACA dataset:

```{r}
system.time({
m = getMACA(AOI = aoi_get(state = "FL"),
model = "CCSM4",
varname = 'pr',
scenario = c('rcp45', 'rcp85'),
startDate = "2080-06-29", endDate = "2080-06-30")
})
```

```{r, fig.width = 15}
ggplot() +
geom_spatraster(data = rast(m)) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white"
) + theme_minimal()
```

Getting multiple models results is also quite simple:

```{r}
models = c("BNU-ESM","CanESM2", "CCSM4")
temp = getMACA(AOI = aoi_get(state = "CO"),
varname = 'tasmin',
model = models,
startDate = "2080-11-29")
temp$air_temperature$mean = terra::app(temp$air_temperature, mean)
names(temp[[1]]) = c(models, "Ensemble Mean")
# Plot
ggplot() +
geom_spatraster(data = temp$air_temperature) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white"
) + theme_minimal()
```
If you don't know your models, you can always grab a random set by specifying a number:

```{r, fig.width= 15}
random = getMACA(aoi_get(state = "MI"),
model = 3,
varname = "pr",
startDate = "2050-10-29")
# Plot
ggplot() +
geom_spatraster(data = random$precipitation) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white"
) +
theme_minimal()
```

# Global Datasets

Not all datasets are USA focused either. TerraClimate offers global, monthly data up to the current year for many variables, and CHIRPS provides daily rainfall data:

```{r, fig.width = 15}
kenya = aoi_get(country = "Kenya")
tc = getTerraClim(kenya, varname = "pet", startDate = "2018-01-01")
chirps = getCHIRPS(kenya, startDate = "2018-01-01", endDate = "2018-01-04" )
library(patchwork)
ggplot() +
geom_spatraster(data = tc$pet) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white"
) +
geom_sf(data = kenya, fill = NA, lwd = 2, col = "black") +
theme_minimal() +
ggplot() +
geom_spatraster(data = chirps$precip) +
facet_wrap(~lyr) +
scale_fill_whitebox_c(
palette = "muted",
na.value = "white"
) +
geom_sf(data = kenya, fill = NA, lwd = 2, col = "black") +
theme_minimal()
```

# Point Based Data

Finally, data gathering is not limited to areal extents and can be retrieved as a time series at locations.

```{r}
AOI =
ts = getGridMET(geocode('Fort Collins', pt = TRUE),
varname = c("pr", 'srad'),
startDate = "2021-01-01",
endDate = "2021-12-31")
```

```{r}
ggplot(data = ts, aes(x = date, y = srad)) +
geom_line() +
stat_smooth(col = "red") +
theme_linedraw() +
labs(title = "Solar Radiation: Fort Collins 2019",
x = "Date", y = "Solar Radiation") +
ggplot(data = ts, aes(x = date, y = cumsum(pr))) +
geom_line(color = "blue", lwd = 1) +
theme_linedraw() +
labs(title = "Cumulative Rainfall: Fort Collins 2019",
x = "Date", y = "Rainfall")
```

# Point Based Ensemble

```{r}
future = getMACA(geocode("Fort Collins", pt = TRUE),
model = 5, varname = "tasmax",
startDate = "2050-01-01", endDate = "2050-01-31")
future_long = pivot_longer(future, -date)
ggplot(data = future_long, aes(x = date, y = value, col = name)) +
geom_line() +
theme_linedraw() +
scale_color_brewer(palette = "Dark2") +
labs(title = "Fort Collins Temperture: January, 2050",
x = "Date",
y = "Degree K",
color = "Model")
```


# Multi site extraction

Extracting data for a set of points is an interesting challenge. It turns it is much more efficient to grab the underlying raster stack and then extract time series as opposed to iterating over the locations:

1. Starting with a set of locations in Colorado:

```{r}
f = system.file("co/co_cities.rds", package = "climateR")
cities = readRDS(f)
```

2. `climateR` will grab the SpatRaster underlying the bounding area of the points
d = getGridMET(AOI,
varname = "pr",
startDate = "1991-10-29",
endDate = "1991-11-06")
```{r}
sites_stack = getTerraClim(AOI = cities,
varname = "tmax",
startDate = "2018-01-01",
endDate = "2018-12-31")
plot(d$precipitation_amount)
```

```{r}
{
plot(sites_stack$tmax[[1]])
plot(vect(cities), add = TRUE, pch = 16, cex = .5)
}
```

3. Use `extract_sites` to extract the times series from these locations. The `id` parameter is the unique identifier from the site data with which to names the resulting columns.

```{r}
sites_wide = extract_sites(r = sites_stack, pts = cities, id = "NAME")
sites_wide[[1]][1:5, 1:5]
```

To make the data 'tidy' simply pivot on the `date` column:

```{r,echo = FALSE}
tmax = tidyr::pivot_longer(sites_wide[[1]], -date)
head(tmax)
ggplot(data = tmax, aes(x = date, y = value, color = name, group = name)) +
scale_color_viridis_d() +
geom_line() +
theme_linedraw() +
theme(legend.position = "none")
```

Loading

0 comments on commit 9aedaf3

Please sign in to comment.