Random fields with physical barriers

class: center, middle, inverse, title-slide

# Random fields with physical barriers
## NOAA PSAW Seminar Series
### March 9, 2022

---

# Adding a barrier to a random field?

Why? Spatial correlation affected by:
  * Coastlines / Islands / lakes ([Bakka 2016](https://arxiv.org/abs/1608.03787))

.center[
<img src="images/bakka_2016.png" width="450px" height = "400px"/>
]

---

# Example of barrier effects on a random field

.center[
<img src="images/barrier.png" width="390px"/>
]

---

# Adding a barrier within INLA

Specify `boundary` attribute to `inla.mesh.2d`

```r
mesh <- inla.mesh.2d(
  boundary = poly,
  ...
)
```

poly created as an `sp` object  
* [INLA book](https://becarioprecario.bitbucket.io/spde-gitbook/ch-nonstationarity.html#ch:barrier)  
* [Haakon Bakka book](https://haakonbakkagit.github.io/btopic128.html)

---

# Example: trawl surveys from Salish Sea

WDFW has extensive survey data on US side of border
  * described and analyzed in [Essington et al. (2021)](https://www.int-res.com/abstracts/meps/v657/p173-189/)

We'll use a single snapshot (2011)
  * focusing on presence-absence of Dungeness crab

---

# Crab data

Presence/absence of Dungeness crabs:

.small[
<img src="11-barrier-models_files/figure-html/load-crabs-1.png" width="700px" style="display: block; margin: auto;" />
]
---

# Get coastline data for Puget Sound

.small[

```r
remotes::install_github("ropensci/rnaturalearthhires")
```

```r
map_data <- rnaturalearth::ne_countries(
  scale = "large",
  returnclass = "sf"
)

puso <- suppressWarnings(suppressMessages(
  sf::st_crop(
    map_data,
    c(xmin = -125, ymin = 46.8, xmax = -122, ymax = 49)
  )
))
crs_utm10 <- 3157 # Pick a projection, here UTM10

# 'WGS84'; necessary on some installs:
sf::st_crs(puso) <- 4326
puso <- sf::st_transform(puso, crs_utm10)
```
]

---

# Project and plot

.small[

```r
survey <- crabs %>%
  sf::st_as_sf(crs = 4326, coords = c("lon", "lat")) %>%
  sf::st_transform(3157) # utm zone 10
ggplot(puso) + # coastline and data:
  geom_sf() +
  geom_sf(data = survey, size = 0.5)
```

<img src="11-barrier-models_files/figure-html/proj-plot-1.png" width="700px" style="display: block; margin: auto;" />
]

---

# Make the mesh

.small[

```r
# Extract the coordinates:
surv_utm_coords <- st_coordinates(survey)

# Scale coordinates to km so the range parameter
# is on a reasonable scale for estimation:
crabs$Xkm <- surv_utm_coords[, 1] / 1000
crabs$Ykm <- surv_utm_coords[, 2] / 1000

mesh <- make_mesh(crabs,
  xy_cols = c("Xkm", "Ykm"),
  cutoff = 7
)
```
]

---

# Plotting the initial mesh

---

# Adding barrier attributes to the mesh

Use `sdmTMB::add_barrier_mesh()` to add on barrier component.

Choose a fraction of the range across physical barriers. Values of 0.1 or 0.2 seem to work well.

```r
barrier_mesh <- add_barrier_mesh(
  spde_obj = mesh,
  barrier_sf = puso,
* range_fraction = 0.1,
  proj_scaling = 1000, # data km but projection m
  plot = FALSE
)
```

---

# Mesh knots/vertices over land (green) vs. water (blue)

---

# Fitting the model

* Fit model as before with a barrier

```r
fit_barrier <- sdmTMB(crab ~ 1,
  data = crabs,
  mesh = barrier_mesh,
  family = binomial(link = "logit")
)
```

---

# Compare spatial parameters

* Without barrier

```
#>      term  estimate std.error  conf.low  conf.high
#> 1   range 31.465826        NA 0.7904659 1252.55011
#> 3 sigma_O  3.574568        NA 0.5005919   25.52485
```

* With barrier

```
#>      term  estimate std.error conf.low conf.high
#> 1   range 20.023808        NA 5.007244  80.07457
#> 3 sigma_O  9.928765        NA 2.360195  41.76789
```

---

# Conclusions

* Accounting for the barrier changes inference about spatial process

In this example:

* Without: less complexity in spatial field (large range, small variance)

* With: larger spatial variation, smaller better estimated range

---

# Advice

* Physical barriers aren't an issue for some coastlines

* Don't choose a `range_fraction` that is much lower than 0.1

* Accounting for spatial correlation barriers can affect local spatial predictions

* Accounting for spatial correlation barriers can affect estimates of range and variance

* But, if the goal is to sum up biomass or abundance across space (e.g., index standardization), the impact is usually subtle