Lecture 4

Introduction to Spatial Statistics

Sara Martino
sara.martino@ntnu.no

Dept. of Mathematical Science, NTNU

Janine Illian
Janine.Illian@glasgow.ac.uk

University of Glasgow

Jafet Belmont
Jafet.Belmont@glasgow.ac.uk

University of Glasgow

Spatial data

Spatial modelling: why necessary?

  • Many natural processes take place in space

  • Spatial data are data which have any form of geographical information attached to them.

  • The emergence of modelling frameworks for spatial data in areas such as epidemiology, ecology and environmental sciences has been facilitated by the rise of new technologies.

Spatial modelling: why necessary?

Large amounts of data collected in space:

increased resolution \(\rightarrow\) large, complex datasets

\(\rightarrow\) complex spatial models required

challenges:

  • often inaccessible to practitioners (literature aimed at statisticians)
  • methodology not always linked to applications
  • models may be too simple to reflect real-life data
  • difficult to apply without expertise in spatial stats and computational statistics

We will see that inlabru can help…

Spatial modelling: why necessary?

Tobler’s first law of geography states that:

Everything is related to everything else, but near things are more related than distant things”

Standard models assume independent observations

  • Spatial data are often not independent

  • Two nearby observations are similar \(\rightarrow\) not providing independent information

\(\rightarrow\) spuriously tight confidence intervals

\(\rightarrow\) wrong inference and wrong conclusions

Spatial models include special components to explicitly model dependence.

Spatial dependence

  • Inference and prediction are essential in statistics, but spatial dependency can complicate statistical analysis.

  • Understanding how spatial dependence structures impact on our analyses is crucial.

  • Spatial dependence may be caused by:

    1. Endogenous processes inherent to the system
      (e.g., localized dispersal leading to organism clustering or social and grouping behaviors)
    2. Exogenous factors such as spatially dependent environmental gradients driving the process of interest
    3. Modelling mis-specifications, including the omission of key covariates or incorrect functional assumptions

Spatial dependence

Dependence relationships need to represented in the model –- which observation is related to which other observations?

  • usually represented by a matrix
  • some matrix operations (e.g. inversion) are computationally expensive
  • in the past: often MCMC ⇝ very slow
  • INLA (and hence inlabru) much faster…

Spatial modelling

The overall goal of any piece of spatial analysis is to understand the spatial patterns in the data. This could involve:

  • Estimating differences in mean, variance, or some other summary statistic across space.

  • Predicting the value at some unobserved location.

  • Identifying hotspots with high (or low) values compared to the rest of the region.

aims also vary with different types of spatial data…

Spatial data structures

Types of spatial data

Discrete space:

  • Data on a spatial grid (areal data)

Continuous space:

  • Geostatistical (geo-referenced) data

  • Spatial point pattern data

The components of the models we will cover are used to reflect spatial dependence structures in discrete and continuous space.

Spatial data structures

We can distinguish three types of spatial data structures

Areal data

Map of bird conservation regions (BCRs) showing the proportion of bird species within each region showing a declining trend

Geostatistical data

Scotland river temperature monitoring network

Point pattern data

Occurrence records of four ungulate species in Tibet,

Spatial data structures

We can distinguish three types of spatial data structures

Areal data

In areal data our measurements are summarised across a set of discrete, non-overlapping spatial units.

Map of bird conservation regions (BCRs) showing the proportion of bird species within each region showing a declining trend

Geostatistical data

Scotland river temperature monitoring network

Point pattern data

Occurrence records of four ungulate species in Tibet

Spatial data structures

We can distinguish three types of spatial data structures

Areal data

Map of bird conservation regions (BCRs) showing the proportion of bird species within each region showing a declining trend

Geostatistical data

In geostatistical data, measurements of a continuous process are taken at a set of fixed locations.

Scotland river temperature monitoring network

Point pattern data

Occurrence records of four ungulate species in Tibet,

Spatial data structures

We can distinguish three types of spatial data structures

Areal data

Map of bird conservation regions (BCRs) showing the proportion of bird species within each region showing a declining trend

Geostatistical data

Scotland river temperature monitoring network

Point pattern data

In point pattern data we record the locations where events occur (e.g. trees in a forest, earthquakes) and the coordinates of such occurrences are our data.

Occurrence records of four ungulate species in Tibet,

Discrete space: areal data

Many public health studies use data aggregated over groups rather than data on individuals - often this is for privacy reasons, but it may also be for convenience.

Respiratory hospitalisations in Glasgow

  • areal data on respiratory hospitalisations
  • 271 Intermediate Zones (IZ) in Greater Glasgow and Clyde health board
  • covariates concern air pollution concentration (PM\(_{10}\)) and socio-economic deprivation in Scotland

Observed response(s): Measurement associated with or areal unit

here: spatial structure represented by Gaussian Markov random field (GMRF)

Figure 1: Greater Glasgow and Clyde health board represented by 271 Intermediate Zones

the data we will use

In this example we model the numbers of hospitalisations due to respiratory disease in the Greater Glasgow and Clyde health board

In epidemiology, disease risk is assessed using Standardized Mortality Ratios (SMR):

\[ SMR_i = \dfrac{Y_i}{E_i} \]

  • A value \(SMR > 1\) indicates a high risk area.
  • A value \(SMR<1\) suggests a low risk area.
  • SMRs may be misleading in counties with small populations.
  • model-based approaches enable us to incorporate covariates and borrow information from neighbouring counties to improve local estimates

Continuous space: geostatistical data

  • We are interested in a phenomenon that is continuous in space
  • The process of interest is measured at a finite set of locations
  • For example, air pollution influencing hospitalisation cases varies continuously in space and does not recognise areal unit boundaries.
  • Our goal is to estimate the value of our variable across the entire space and/or to model the relationship with covariates.

Observed response(s): measurements at given locations in continuous space

here: spatial structure represented by Gaussian random field; approximated by a continuously indexed Gauss Markov random field

the data we will use

Data on the Pacific Cod (Gadus macrocephalus) from a trawl survey in Queen Charlotte Sound.

  • Presence/Absence of Pacific cod in the area swept for a given survey in 2003 as well as depth covariate information.

Continuous space: spatial point patterns

  • Locations of objects/events in space (typically 2D)
  • Examples: tree locations, animal groups, earthquakes

Observed response(s): x,y coordinates (sometimes also additional measurements; “marks”) modelled by a random variable, a spatial point process characterised by intensity \(\lambda(s) \in \mathbb{R}^2\)

here: intensity is a Gaussian random field; approximated by a continuously indexed Gauss Markov random field

the data we will use

locations of forest fires in the Castilla-La Mancha region of Spain between 1998 and 2007

citizen science records of the ringlet butterfly in Scotland’s Cairngorms National Park

Spatial Data and inlabru

inlabru provides support for sf (Pebesma, 2018) and terra (Hijmans, 2024) spatial data structures.

  • sf and terra objects can be passed as data and covariates directly to bru_obs() and bru() functions.

In this course we will

  • Understand how eval_spatial() methods are used to extract information from spatial data objects.

  • Define model components based on spatial information provided in the geometry column of the sf data objects.

  • Learn about the bru_mapper() mapper system for specifying SPDE random effect.

  • Use Spatial objects for prediction using predict() and generate().

  • Fit log-Gaussian Cox process using (i) sf points object to describe the locations of the observed points and (ii) sf polygon object to define the observation window

First, we will cover some of the basics of working sf and terra objects.

A very short overview of sf and terra

Why these two packages?

Modern, efficient framework for working with spatial data.

Our spatial analysis begins with how the data are stored:

  • sf (Simple Features ) package → vector data. Discrete features with exact coordinates: points (a survey location), lines (a transect), polygons (a health board, a national park).

  • It can also be used to read vector data stored as a shapefiles.

  • terra → raster data. A continuous surface(s) split into a regular grid, each holding a value: elevation, sea depth, temperature.

sf: vector data as a data frame

The key idea: an sf object is a data frame — with one extra geometry column.

Three functions do most of the work:

  • st_read() — read a shapefile (.shp) from disk
  • st_as_sf() — turn a data frame of coordinates into spatial points
  • st_transform() — reproject to a different coordinate system
my_sf <- st_as_sf(my_df, coords = c("lon", "lat"), crs = 4326)

Plot it by adding geom_sf() layer to a ggplot() object.

sf: vector data as a data frame

Most Tidyverse/dplyr verbs (e.g., filter, mutate,summarise, %>%) work with sf object, providing a seamless bridge between data manipulation and statistical analysis.

Coordinate Reference System

A Coordinate Reference System tells R what your coordinates mean — degrees of latitude/longitude, or metres on a projected grid. E.g.,

  • EPSG:4326 — raw longitude/latitude (WGS84), in degrees
  • UTM zones (e.g. EPSG:32609, 27700) — projected, in metres, better for measuring distance and area
  1. Layers you overlay must share the same CRS — reproject using st_transform()
  2. Switching units from metres to kilometres makes areas and distances easier to read and interpret.

Important

When your spatial data carry a CRS, inlabru uses it directly — so make sure all layers in your data share a consistent and appropriate CRS before you fit. Mismatched coordinate systems are one of the most common sources of error.

terra: raster data

Rasters are a convenient way to represents spatially continuous phenomena by dividing a region into a grid of equally-sized cells, each storing a value for the variable of interest

The terra package is a modern and powerful tool for efficiently working with raster data.

Some useful functions are:

  • rast() — read a raster file (.tif or .tiff) or build one from an x-y-z data frame
  • crs() — assign or match a coordinate system
  • crop() — clip a raster to a region of interest
  • ext() - get the spatial extension
  • extend()- Enlarge the spatial extent of a raster
raster_data <- rast(qcs_grid, type = "xyz")
crs(raster_data) <- crs(sf_object)   # match the vector layer

To plot a raster in ggplot, the tidyterra package adds geom_spatraster() — so raster and sf layers stack on the same figure.

Take Home message

  • all spatial models discussed are special cases of Latent Gaussian models

  • different spatial terms are needed for different spatial data structures

  • the SPDE approach unifies approximation as all data structures are spatially referenced in continuous space

  • inlabru fits these models efficiently

In the next practical you will explore tools for visualization and wrangling spatial data objects.