Quetzal-CoaTL
The Coalescence Template Library
Loading...
Searching...
No Matches
Autoregressive processes and time series

Background

While local parameters describe a static relationship to the environment, local processes on a spatial graph are dynamic in nature. They encompass the dynamic phenomena occurring within localized regions defined by the vertices of the graph, including interactions, transformations, and behaviors unfolding at a granular level. Understanding these processes is vital for comprehending how spatial relationships evolve, how information or resources propagate through the graph, and how localized phenomena contribute to broader spatial patterns and their temporal trends.

In the realm of ecology and evolution, these local processes encompass the intricate dynamics unfolding within specific ecological niches or habitats delineated by the vertices of the graph. These processes encapsulate ecological interactions, evolutionary adaptations, and population dynamics occurring at localized scales, shaping the distribution, abundance, and genetic diversity of species within the spatial landscape.

Formalizing these local processes is crucial for understanding ecological community assembly, species responses to environmental changes, and the emergence of biodiversity patterns in spatially diverse landscapes. By investigating local ecological and evolutionary dynamics, researchers can uncover the mechanisms underlying species interactions, habitat fragmentation, and speciation events, thus deepening our understanding of ecosystem functioning and resilience in the face of environmental perturbations.

Local population growth models, integral to ecological and evolutionary studies, allow researchers to simulate and comprehend population dynamics within specific geographic areas or habitats. These models consider environmental conditions, species interactions, and demographic processes to forecast changes in population size and composition over time.

This page provides examples of how to combine Quetzal components to construct the local simulation process corresponding to your hypotheses.

Autoregressive Models

A process in which the value at time \( t \) depends on its value at time \( t−1 \) (like population growth over time) is commonly referred to as an autoregressive process or a first-order autoregressive process ( \(AR(1)\)). These are widely used in time series analysis and modeling to capture temporal dependencies and patterns in data.

In mathematical terms, it can be represented as:

\[ X_t =​ \phi X_{t-1} + \epsilon_t \]

where:

  • \( X_t \) is the value of the process at time \( t \)
  • \( X_{t-1} \) is the value of the process at time \( t-1 \)
  • \( \phi \) is the autoregressive parameter,
  • \( \epsilon_t \) is a random error term.

Users can represent first-order (or even p-order) autoregressive processes using Quetzal components. Although flexible, the general approach will be to represent time series as std::vector structures embedded in each vertex or edges they relate to, where a quantity \(q\) a time \(t\) at vertex \(x\) is described by graph[x][q][t].


Deterministic Models

Deterministic models describe population growth based on precise mathematical equations that govern population dynamics. These models typically assume constant parameters and predictable relationships between population variables, such as birth rates, death rates, and carrying capacity. While deterministic models provide valuable insights into long-term population trends, they may oversimplify real-world complexities and fail to account for stochastic fluctuations or environmental variability.

In this example, we store a time-series vector at each vertex of the graph to monitor the population dynamics. After initializing one vertex with a non-zero value, we iterate through time to update the series. It's important to note that, since we defined no migration along the graph's edges, all time-series except one remain at zero.

Input

1#include "quetzal/geography.hpp"
2
3namespace geo = quetzal::geography;
4
5int main()
6{
7 auto file = std::filesystem::current_path() / "data/bio1.tif";
8
9 // The raster have 10 bands that we will assign to 2001 ... 2010.
10 std::vector<int> times(10);
11 std::iota(times.begin(), times.end(), 2001);
12
13 // Landscape construction injecting times information
14 auto land = geo::raster<>::from_file(file, times);
15
16 // Graph construction with vertex and edge info
17 using vertex_info = std::vector<double>;
18 vertex_info N(times.size(), 0);
19
22 graph[0][0] = 1;
23
24 // Part of what happens in a demographic simulator:
25 for( size_t t = 1; t < times.size(); ++t){
26 for( auto x : graph.vertices() ){
27 graph[x][t] = 2 * graph[x][t-1];
28 }
29 }
30
31 // Post treatment
32 std::cout << "Time series:\n";
33 for(auto x : graph.vertices()){
34 std::cout << "at vertex " << x << " : ";
35 for(auto value : graph[x] ){
36 std::cout << value << " ";
37 }
38 std::cout << "\n";
39 }
40}
Individuals can not escape the landscape's borders.
Definition bound_policy.hpp:20
Discrete spatio-temporal variations of an environmental variable.
Definition raster.hpp:38
Geospatial data formatting and processing.
Definition geography.hpp:17
boost::no_property no_property
Represents no information carried by vertices or edges of a graph.
Definition concepts.hpp:71
auto from_grid(SpatialGrid const &grid, VertexProperty const &v, EdgeProperty const &e, Vicinity const &vicinity, Directionality dir, Policy const &bounding_policy)
Spatial graph construction method.
Definition from_grid.hpp:36
boost::undirectedS isotropy
Property of a process independent of the direction of movement.
Definition directionality.hpp:18
Definition geography_dispersal_kernel_4.cpp:13
Definition vicinity.hpp:53
Definition coalescence_binary_tree_2.cpp:5

Output

1Time series:
2at vertex 0 : 1 2 4 8 16 32 64 128 256 512
3at vertex 1 : 0 0 0 0 0 0 0 0 0 0
4at vertex 2 : 0 0 0 0 0 0 0 0 0 0
5at vertex 3 : 0 0 0 0 0 0 0 0 0 0
6at vertex 4 : 0 0 0 0 0 0 0 0 0 0
7at vertex 5 : 0 0 0 0 0 0 0 0 0 0
8at vertex 6 : 0 0 0 0 0 0 0 0 0 0
9at vertex 7 : 0 0 0 0 0 0 0 0 0 0
10at vertex 8 : 0 0 0 0 0 0 0 0 0 0

Stochastic Models

Stochastic models introduce randomness and variability into population growth processes, acknowledging the inherent uncertainty in ecological systems. These models incorporate probabilistic elements, such as random fluctuations in birth and death rates, environmental variability, and demographic stochasticity. Stochastic models are particularly useful for capturing the effects of environmental variability and demographic stochasticity on population dynamics, allowing researchers to assess the likelihood of rare events and population extinctions.

We build up on the previous example by adding some stochasticity in the process. We store a time-series vector at each vertex of the graph to monitor the population dynamics. After initializing some vertices with a non-zero value, we iterate through time to update the series. It's important to note that, since we defined no migration along the graph's edges, some time-series remain at zero.

Input

1#include "quetzal/geography.hpp"
2#include <random>
3#include <cmath>
4
5namespace geo = quetzal::geography;
6
7int main()
8{
9 auto file = std::filesystem::current_path() / "data/bio1.tif";
10
11 // The raster have 10 bands that we will assign to 2001 ... 2010.
12 std::vector<int> times(10);
13 std::iota(times.begin(), times.end(), 2001);
14
15 // Landscape construction injecting times information
16 auto land = geo::raster<>::from_file(file, times);
17
18 // Graph construction with vertex and edge info
19 using vertex_info = std::vector<double>;
20 vertex_info N(times.size(), 0);
21
24
25 // Let's compare some time series
26 graph[0][0] = 1;
27 graph[1][0] = 1;
28 graph[2][0] = 1;
29
30 // Declare a random number generator
31 std::mt19937 rng{std::random_device{}()};
32
33 // Part of what happens in a demographic simulator:
34 for( size_t t = 1; t < times.size(); ++t){
35 for( auto x : graph.vertices() ){
36 auto alpha = std::uniform_real_distribution<double>(1, 2)(rng);
37 graph[x][t] = alpha * graph[x][t-1];
38 }
39 }
40
41 // Post treatment
42 std::cout << "Time series:\n";
43 for(auto x : graph.vertices()){
44 std::cout << "at vertex " << x << " : ";
45 for(auto value : graph[x] ){
46 std::cout << static_cast<int>(value) << " ";
47 }
48 std::cout << "\n";
49 }
50}

Output

1Time series:
2at vertex 0 : 1 1 1 2 3 5 9 14 23 40
3at vertex 1 : 1 1 2 4 8 12 19 34 59 63
4at vertex 2 : 1 1 1 1 2 4 6 9 11 21
5at vertex 3 : 0 0 0 0 0 0 0 0 0 0
6at vertex 4 : 0 0 0 0 0 0 0 0 0 0
7at vertex 5 : 0 0 0 0 0 0 0 0 0 0
8at vertex 6 : 0 0 0 0 0 0 0 0 0 0
9at vertex 7 : 0 0 0 0 0 0 0 0 0 0
10at vertex 8 : 0 0 0 0 0 0 0 0 0 0

Environmental Niche Models

Environmental niche models focus on how environmental factors influence the distribution and abundance of species within specific habitats. These models integrate ecological niche theory and environmental data to predict species' occurrence probabilities and population densities across spatial gradients. By quantifying the relationships between environmental variables (e.g., temperature, precipitation, habitat quality) and species' ecological requirements, environmental niche models help elucidate how environmental factors shape local population growth and species distributions.

In this example, we store a time-series vector at each vertex of the graph to monitor the population dynamics. After initializing some vertices with a non-zero value, we iterate through time to update the series conditionally to the local contemporary environmental data. It's important to note that, since we defined no migration along the graph's edges, some time-series remain at zero.

Since the data structures embedded in the vertices (and edges) of the graph are arbitrary, we could also store the local \(K\) and \(r\) time series alongside the population size. This can be achieved by defining a simple vertex_info class that contains three std::vector members. We will apply this idea in the following example for simulating a three-species competition model.

Input

1#include "quetzal/geography.hpp"
2#include <random>
3#include <cmath>
4
5namespace geo = quetzal::geography;
6
7int main()
8{
9 auto file = std::filesystem::current_path() / "data/bio1.tif";
10
11 // The raster have 10 bands that we will assign to 2001 ... 2010.
12 std::vector<int> times(10);
13 std::iota(times.begin(), times.end(), 2001);
14
15 // Landscape construction injecting times information
16 auto land = geo::raster<>::from_file(file, times);
17
18 // Graph construction with vertex and edge info
19 using vertex_info = std::vector<double>;
20 vertex_info N(times.size(), 0);
21
24
25 // Let's compare some time series
26 graph[0][0] = 1;
27 graph[1][0] = 1;
28 graph[2][0] = 1;
29
30 // Declare a random number generator
31 std::mt19937 rng{std::random_device{}()};
32
33 // Part of what could happen in a demographic simulator:
34 for( size_t t = 1; t < times.size(); ++t){
35 for( auto x : graph.vertices() ){
36 auto r = std::uniform_real_distribution<double>(1, 5)(rng);
37 auto K = 10. * land.at(x, t-1).value_or(10.);
38 auto P = graph[x][t-1];
39 graph[x][t] = std::max(P + r * P * (1 - P/K), 0.);
40 }
41 }
42
43 // Post treatment
44 std::cout << "Time series:\n";
45 for(auto x : graph.vertices()){
46 std::cout << "at vertex " << x << " : ";
47 for(auto value : graph[x] ){
48 std::cout << static_cast<int>(value) << " ";
49 }
50 std::cout << "\n";
51 }
52}

Output

1Time series:
2at vertex 0 : 1 3 20 99 432 1448 358 788 1690 0
3at vertex 1 : 1 2 9 51 272 1257 616 1616 774 1772
4at vertex 2 : 1 4 12 64 323 935 1088 969 1135 942
5at vertex 3 : 0 0 0 0 0 0 0 0 0 0
6at vertex 4 : 0 0 0 0 0 0 0 0 0 0
7at vertex 5 : 0 0 0 0 0 0 0 0 0 0
8at vertex 6 : 0 0 0 0 0 0 0 0 0 0
9at vertex 7 : 0 0 0 0 0 0 0 0 0 0
10at vertex 8 : 0 0 0 0 0 0 0 0 0 0

Species Competition Models

Species competition models explore how interspecific interactions, such as competition for resources or interactions with predators, influence population dynamics and species coexistence within local communities. These models typically incorporate competition coefficients, representing the strength of interspecific interactions, into population equations to simulate the effects of species interactions on population growth rates and carrying capacities. Species competition models provide insights into the mechanisms driving species coexistence, competitive exclusion, and community dynamics within ecological communities.

A discrete-time model for a three-species predator-prey competition can be complex, but here's a basic version using a system of difference equations. Let's denote the populations of the three species at time \(t\) as \(𝑃_1(t)\), \(P_2(t)\) and \(P_3(t)\). Assume that \(P_1\) is a prey species, \(P_2\) is a predator of \(P_1\) and \(P_3\) is a competing species that preys on \(P_1\) and competes with \(P_2\).

\[ P_1(t+1) = P_1(t) + r_1 P_1(t) \left(1 - \frac{P_1(t)}{K_1}\right) - a_{12} P_1(t) P_2(t) - a_{13} P_1(t) P_3(t) \]

\[ P_2(t+1) = P_2(t) + r_2 P_2(t) \left(1 - \frac{P_2(t)}{K_2}\right) + b_{21} P_1(t) P_2(t) - c_{23} P_2(t) P_3(t) \]

\[ P_3(t+1) = P_3(t) + r_3 P_3(t) \left(1 - \frac{P_3(t)}{K_3}\right) + b_{31} P_1(t) P_3(t) - c_{32} P_2(t) P_3(t) \]

Where:

  • \( r_1, r_2, r_3 \) are the intrinsic growth rates of species 1, 2, and 3, respectively.
  • \( K_1, K_2, K_3 \) are the carrying capacities for species 1, 2, and 3, respectively.
  • \( a_{12} \) is the predation rate of species 2 on species 1.
  • \( a_{13} \) is the predation rate of species 3 on species 1.
  • \( b_{21} \) is the benefit species 2 gets from preying on species 1.
  • \( c_{23} \) is the competition rate between species 2 and 3.
  • \( b_{31} \) is the benefit species 3 gets from preying on species 1.
  • \( c_{32} \) is the competition rate between species 3 and 2.

Input

1#include "quetzal/geography.hpp"
2#include <random>
3#include <cmath>
4
5namespace geo = quetzal::geography;
6
7// Time series of the 3 species population sizes to embed along the graph
8struct vertex_info {
9 using time_series_type = std::vector<double>;
10 time_series_type P1, P2, P3;
11 // Constructor to initialize the vectors with a specified length and fill them with zeros
12 vertex_info(int length): P1(length, 0), P2(length, 0), P3(length, 0){}
13 // Needs to be default constructible
14 vertex_info() = default;
15};
16
17// Overload the << operator for the vertex_info struct
18std::ostream& operator<<(std::ostream& os, const vertex_info& vertex) {
19 auto format = [&](const auto & v, const auto& name){
20 os << "\t" << name << " : ";
21 for(const auto& val : v)
22 os << val << " ";
23 os << "\n";
24 };
25 format(vertex.P1, "P1");
26 format(vertex.P2, "P2");
27 format(vertex.P3, "P3");
28 return os;
29}
30
31int main()
32{
33 auto file = std::filesystem::current_path() / "data/bio1.tif";
34
35 // The raster have 10 bands that we will assign to 2001 ... 2010.
36 std::vector<int> times(10);
37 std::iota(times.begin(), times.end(), 2001);
38
39 // Landscape construction injecting times information
40 auto land = geo::raster<>::from_file(file, times);
41
42 // Graph construction with vertex and edge info
44 auto graph = geo::from_grid(land, vertex_info(times.size()), edge_info(), geo::connect_fully(), geo::isotropy(), geo::mirror());
45
46 // Let's initialize population sizes for a single vertex for demonstration
47 graph[0].P1[0] = 200;
48 graph[0].P2[0] = 1;
49 graph[0].P3[0] = 1;
50
51 // Declare a random number generator
52 std::mt19937 rng{std::random_device{}()};
53
54 auto random_double = [&](double min, double max){ return std::uniform_real_distribution<double>(min, max)(rng);};
55
56 double r1 = random_double(0.0, 1.0); // Intrinsic growth rate for species 1
57 double r2 = random_double(0.0, 1.0); // Intrinsic growth rate for species 2
58 double r3 = random_double(0.0, 1.0); // Intrinsic growth rate for species 3
59
60 double K1 = random_double(50.0, 100.0); // Carrying capacity for species 1
61 double K2 = random_double(50.0, 100.0); // Carrying capacity for species 2
62 double K3 = random_double(50.0, 100.0); // Carrying capacity for species 3
63
64 double a12 = random_double(0.0, 0.1); // Predation rate of species 2 on species 1
65 double a13 = random_double(0.0, 0.1); // Predation rate of species 3 on species 1
66
67 double b21 = random_double(0.0, 0.1); // Benefit species 2 gets from preying on species 1
68 double c23 = random_double(0.0, 0.1); // Competition rate between species 2 and 3
69
70 double b31 = random_double(0.0, 0.1); // Benefit species 3 gets from preying on species 1
71 double c32 = random_double(0.0, 0.1); // Competition rate between species 3 and 2
72
73 // Part of what could happen in a demographic simulator:
74 for( size_t t = 1; t < times.size(); ++t){
75 for( auto x : graph.vertices() ){
76 auto P1 = graph[x].P1[t-1];
77 auto P2 = graph[x].P2[t-1];
78 auto P3 = graph[x].P3[t-1];
79 graph[x].P1[t] = std::max(P1 + r1 * P1 * (1 - P1 / K1) - a12 * P1 * P2 - a13 * P1 * P3, 0.);
80 graph[x].P2[t] = std::max(P2 + r2 * P2 * (1 - P2 / K2) + b21 * P1 * P2 - c23 * P2 * P3, 0.);
81 graph[x].P3[t] = std::max(P3 + r3 * P3 * (1 - P3 / K3) + b31 * P1 * P3 - c32 * P2 * P3, 0.);
82 }
83 }
84
85 // Post treatment
86 std::cout << "Time series:\n";
87 for(auto x : graph.vertices()){
88 std::cout << "at vertex " << x << " : \n";
89 std::cout << graph[x] << std::endl;
90 }
91}

Output

1Time series:
2at vertex 0 :
3 P1 : 200 0 0 0 0 0 0 0 0 0
4 P2 : 1 9.64755 12.668 15.1625 16.2632 15.6243 13.4454 9.96643 5.48145 1.08938
5 P3 : 1 2.07459 3.12903 4.30814 5.44838 6.58124 8.02821 10.3471 14.4745 22.1089
6
7at vertex 1 :
8 P1 : 0 0 0 0 0 0 0 0 0 0
9 P2 : 0 0 0 0 0 0 0 0 0 0
10 P3 : 0 0 0 0 0 0 0 0 0 0
11
12at vertex 2 :
13 P1 : 0 0 0 0 0 0 0 0 0 0
14 P2 : 0 0 0 0 0 0 0 0 0 0
15 P3 : 0 0 0 0 0 0 0 0 0 0
16
17at vertex 3 :
18 P1 : 0 0 0 0 0 0 0 0 0 0
19 P2 : 0 0 0 0 0 0 0 0 0 0
20 P3 : 0 0 0 0 0 0 0 0 0 0
21
22at vertex 4 :
23 P1 : 0 0 0 0 0 0 0 0 0 0
24 P2 : 0 0 0 0 0 0 0 0 0 0
25 P3 : 0 0 0 0 0 0 0 0 0 0
26
27at vertex 5 :
28 P1 : 0 0 0 0 0 0 0 0 0 0
29 P2 : 0 0 0 0 0 0 0 0 0 0
30 P3 : 0 0 0 0 0 0 0 0 0 0
31
32at vertex 6 :
33 P1 : 0 0 0 0 0 0 0 0 0 0
34 P2 : 0 0 0 0 0 0 0 0 0 0
35 P3 : 0 0 0 0 0 0 0 0 0 0
36
37at vertex 7 :
38 P1 : 0 0 0 0 0 0 0 0 0 0
39 P2 : 0 0 0 0 0 0 0 0 0 0
40 P3 : 0 0 0 0 0 0 0 0 0 0
41
42at vertex 8 :
43 P1 : 0 0 0 0 0 0 0 0 0 0
44 P2 : 0 0 0 0 0 0 0 0 0 0
45 P3 : 0 0 0 0 0 0 0 0 0 0
46