Merge pull request #1 from EstebanAce/master

update build system and merge bounded_dispersal

.git_archival.txt

@@ -0,0 +1,4 @@
+node: $Format:%H$
+node-date: $Format:%cI$
+describe-name: $Format:%(describe:tags=true,match=*[0-9]*)$
+ref-names: $Format:%D$

.gitattributes

@@ -0,0 +1,3 @@
+adascape/_version.py export-subst
+
+.git_archival.txt  export-subst

.github/workflows/test.yml

@@ -0,0 +1,37 @@
+name: test build
+
+on:
+  push:
+    branches:
+      - master
+  pull_request:
+    branches:
+      - master
+
+jobs:
+  test:
+    name: Run pytest
+    runs-on: ubuntu-latest
+    steps:
+      - name: Checkout
+        uses: actions/checkout@v4
+      - name: Setup micromamba
+        uses: mamba-org/setup-micromamba@v1
+        with:
+          cache-environment: true
+          cache-downloads: false
+          environment-file: environment.yml
+      - name: Conda info
+        shell: bash -l {0}
+        run: conda info
+      - name: Conda list
+        shell: bash -l {0}
+        run: conda list
+      - name: Install adascape
+        shell: bash -l {0}
+        run: |
+          pip install --no-deps -e .
+      - name: Run tests
+        shell: bash -l {0}
+        run: pytest adascape -v
+

.github/workflows/test_notebooks.yml

@@ -0,0 +1,33 @@
+name: test notebooks
+
+on:
+  push:
+    branches:
+    - master
+  pull_request:
+    branches:
+    - master
+
+jobs:
+  run_all_notebooks:
+    name: Run all notebooks
+    runs-on: ubuntu-latest
+    steps:
+      - name: Checkout
+        uses: actions/checkout@v4
+      - name: Setup micromamba
+        uses: mamba-org/setup-micromamba@v1
+        with:
+          cache-environment: true
+          cache-downloads: false
+          environment-file: environment.yml
+      - name: Conda info
+        shell: bash -l {0}
+        run: conda info
+      - name: Install adascape
+        shell: bash -l {0}
+        run: |
+          pip install --no-deps -e .
+      - name: Execute all notebooks
+        shell: bash -l {0}
+        run: python execute_all_notebooks.py

README.md

@@ -1,17 +1,26 @@
+[![test build](https://github.com/fastscape-lem/adascape/actions/workflows/test.yml/badge.svg?branch=master)](https://github.com/fastscape-lem/adascape/actions)
+[![test notebooks](https://github.com/EstebanAce/adascape/actions/workflows/test_notebooks.yml/badge.svg?branch=master)](https://github.com/EstebanAce/adascape/actions)
+[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.7794374.svg)](https://doi.org/10.5281/zenodo.7794374)
+
+
+
 # AdaScape: Adaptive speciation and landscape evolution model
 
-The "AdaScape" package contains a simple adaptive speciation models written
-in Python that can be easily coupled with a landscape evolution model (FastScape).
+The "AdaScape" package contains a simple adaptive speciation model written
+in Python that is coupled with the landscape evolution model [FastScape](https://fastscape.readthedocs.io/en/latest/).
 
 ## Install
 
-This package depends on Python (3.5 or later is recommended),
+This package depends on Python (3.9 or later is recommended),
 [numpy](http://www.numpy.org/),
-[scipy](https://docs.scipy.org/doc/scipy/reference/) and
-[pandas](https://pandas.pydata.org/).
-
-This package also provides a [fastscape](https://fastscape.readthedocs.io)
-and a [dendropy](https://dendropy.org/) extensions (optional dependencies).
+[scipy](https://docs.scipy.org/doc/scipy/reference/),
+[pandas](https://pandas.pydata.org/),
+[fastscape](https://github.com/fastscape-lem/fastscape) and 
+[orographic precipitation](https://github.com/fastscape-lem/orographic-precipitation) .
+
+This package also provides a [dendropy](https://dendropy.org/) extension and 
+uses [toytree](https://toytree.readthedocs.io/en/latest/index.html) 
+to plot phylogenetic trees (optional dependencies).
 
 To install the package locally, first clone this repository:
 

adascape/base.py

@@ -1,3 +1,4 @@
+import pdb
 import textwrap
 import warnings
 import numpy as np
@@ -368,9 +369,10 @@ def _sample_in_range(self, values_range):
         """
         return self._rng.uniform(values_range[0], values_range[1], self._init_abundance)
 
-    def _mov_within_bounds(self, x, y, sigma):
+    def _mov_within_bounds(self, x, y, sigma, disp_boundary=None):
         """
-        Move and check if the location of individuals are within grid range.
+        Move and check if the location of individuals are within grid range
+        or available are to disperse.
 
         Parameters
         ----------
@@ -380,6 +382,9 @@ def _mov_within_bounds(self, x, y, sigma):
             locations along the y coordinate
         sigma : float
             movement variability
+        disp_boundary: float, array-like
+            dispersal boundaries as an array or list of vertices [[x,y],...]
+            delimiting the area where individuals can disperse
         Returns
         -------
         array-like
@@ -391,6 +396,17 @@ def _mov_within_bounds(self, x, y, sigma):
         delta_bounds_y = self._grid_bounds['y'][:, None] - y
         new_x = self._truncnorm.rvs(*(delta_bounds_x / sigma), loc=x, scale=sigma)
         new_y = self._truncnorm.rvs(*(delta_bounds_y / sigma), loc=y, scale=sigma)
+
+        if disp_boundary is not None:
+            hull = spatial.Delaunay(disp_boundary)
+            inhull = hull.find_simplex(np.column_stack([new_x, new_y])) <= 0
+
+            while any(inhull):
+                delta_bounds_x = self._grid_bounds['x'][:, None] - x[inhull]
+                delta_bounds_y = self._grid_bounds['y'][:, None] - y[inhull]
+                new_x[inhull] = self._truncnorm.rvs(*(delta_bounds_x / sigma), loc=x[inhull], scale=sigma)
+                new_y[inhull] = self._truncnorm.rvs(*(delta_bounds_y / sigma), loc=y[inhull], scale=sigma)
+                inhull = hull.find_simplex(np.column_stack([new_x, new_y])) <= 0
         return new_x, new_y
 
     def _mutate_trait(self, trait, sigma):
@@ -412,20 +428,23 @@ def _mutate_trait(self, trait, sigma):
         mut_trait = self._truncnorm.rvs(a, b, loc=trait, scale=sigma)
         return mut_trait
 
-    def _update_individuals(self, dt):
+    def _update_individuals(self, dt, disp_boundary=None):
         """Require implementation in subclasses."""
         raise NotImplementedError()
 
-    def update_individuals(self, dt):
+    def update_individuals(self, dt, disp_boundary=None):
         """Update individuals' data (generate, mutate, and disperse).
 
         Parameters
         ----------
         dt : float
             Time step duration.
+        disp_boundary: float, array-like
+            dispersal boundaries as an array or list of vertices [[x,y],...]
+            delimiting the area where individuals can disperse
 
         """
-        self._update_individuals(dt)
+        self._update_individuals(dt, disp_boundary)
 
         if not self._params['always_direct_parent']:
             self._set_direct_parent = False
@@ -552,15 +571,10 @@ def _count_neighbors(self, locations, traits):
         n_all = np.array([len(nb) for nb in neighbors])
         return n_all, n_eff
 
-    def evaluate_fitness(self, dt):
+    def evaluate_fitness(self):
         """Evaluate fitness and generate offspring number for group of individuals and
         with environmental conditions both taken at the current time step.
 
-        Parameters
-        ----------
-        dt : float
-            Time step duration.
-
         """
 
         if self.abundance:
@@ -602,13 +616,16 @@ def evaluate_fitness(self, dt):
             'n_eff': n_eff
         })
 
-    def _update_individuals(self, dt):
+    def _update_individuals(self, dt, disp_boundary=None):
         """Update individuals' data (generate, mutate, and disperse).
 
         Parameters
         ----------
         dt : float
             Time step duration.
+        disp_boundary: float, array-like
+            dispersal boundaries as an array or list of vertices [[x,y],...]
+            delimiting the area where individuals can disperse
         """
 
         n_offspring = self._individuals['n_offspring']
@@ -644,7 +661,8 @@ def _update_individuals(self, dt):
             # disperse offspring within grid bounds
             new_x, new_y = self._mov_within_bounds(new_individuals['x'],
                                                    new_individuals['y'],
-                                                   self._params['sigma_d'])
+                                                   self._params['sigma_d'],
+                                                   disp_boundary)
             new_individuals['x'] = new_x
             new_individuals['y'] = new_y
 

adascape/fastscape_ext.py

@@ -47,6 +47,11 @@ class Speciation:
     grid_x = xs.foreign(UniformRectilinearGrid2D, "x")
     grid_y = xs.foreign(UniformRectilinearGrid2D, "y")
 
+    disp_boundary = xs.variable(default=None, description="dispersal boundaries as an xr.DataArray "
+                                                          "with vertices [[x,y],...] of bounded area "
+                                                          "with dimensions p and d",
+                                static=True, dims=[(), ('p', 'd')])
+
     _model = xs.any_object(description="speciation model instance")
     _individuals = xs.any_object(description="speciation model state dictionary")
 
@@ -158,20 +163,19 @@ def initialize(self):
 
         self._model.initialize(init_x_range, init_y_range)
 
-    @xs.runtime(args='step_delta')
-    def run_step(self, dt):
+    def run_step(self):
         # reset individuals "cache"
         self._individuals = None
 
         # maybe update model parameters
         self._model.params.update(self._get_model_params())
 
         self.abundance = self._model.abundance
-        self._model.evaluate_fitness(dt)
+        self._model.evaluate_fitness()
 
     @xs.runtime(args='step_delta')
     def finalize_step(self, dt):
-        self._model.update_individuals(dt)
+        self._model.update_individuals(dt, self.disp_boundary)
 
     @fitness.compute
     def _get_fitness(self):

adascape/tests/test_base.py

@@ -209,7 +209,7 @@ def test_evaluate_fitness(self, model, model_IR12):
         if model == 'IR12':
             model_IR12.initialize()
             init_pop = model_IR12.individuals.copy()
-            model_IR12.evaluate_fitness(1)
+            model_IR12.evaluate_fitness()
             eval_pop = model_IR12.individuals.copy()
             for k in ['fitness', 'n_offspring']:
                 assert k in eval_pop
@@ -226,7 +226,7 @@ def test_update_individuals(self, model_IR12, grid):
 
             model_IR12.initialize()
             init_pop = model_IR12.individuals.copy()
-            model_IR12.evaluate_fitness(1)
+            model_IR12.evaluate_fitness()
             model_IR12.update_individuals(1)
             current_pop = model_IR12.individuals.copy()
 
@@ -239,7 +239,7 @@ def test_update_individuals(self, model_IR12, grid):
             assert _in_bounds(grid[1], current_pop['y'])
 
             # test mutation
-            model_IR12.evaluate_fitness(1)
+            model_IR12.evaluate_fitness()
             model_IR12.update_individuals(1)
             last_pop = model_IR12.individuals.copy()
             idx = np.searchsorted(current_pop['taxon_id'], last_pop['ancestor_id'])-1
@@ -248,7 +248,7 @@ def test_update_individuals(self, model_IR12, grid):
         trait_diff = np.concatenate(trait_diff)
         trait_rms = np.sqrt(np.mean(trait_diff ** 2))
         scaled_sigma_m = model_IR12.params['sigma_m'] * np.sqrt(model_IR12.params['p_m'])
-        assert trait_rms == pytest.approx(scaled_sigma_m, 0.1, 0.02)
+        assert trait_rms == pytest.approx(scaled_sigma_m, 0.1, 0.05)
 
         # test reset fitness data
         for k in ['fitness', 'n_offspring']:
@@ -263,19 +263,19 @@ def test_updade_population_parents(self, grid, params_IR12, trait_funcs, direct_
         model = IR12SpeciationModel(X, Y, init_trait_funcs, opt_trait_funcs, 10, **params_IR12)
         model.initialize()
 
-        model.evaluate_fitness(1)
+        model.evaluate_fitness()
         parents0 = model.to_dataframe(varnames='ancestor_id')
         model.update_individuals(1)
 
-        model.evaluate_fitness(1)
+        model.evaluate_fitness()
         parents1 = model.to_dataframe(varnames='ancestor_id')
         model.update_individuals(1)
 
-        model.evaluate_fitness(1)
+        model.evaluate_fitness()
         parents2 = model.to_dataframe(varnames='ancestor_id')
         model.update_individuals(1)
 
-        model.evaluate_fitness(1)
+        model.evaluate_fitness()
         parents3 = model.to_dataframe(varnames='ancestor_id')
         model.update_individuals(1)
 
@@ -311,24 +311,24 @@ def get_pop_subset():
 
         if on_extinction == 'raise':
             with pytest.raises(RuntimeError, match="no offspring"):
-                initialized_model_IR12.evaluate_fitness(1)
+                initialized_model_IR12.evaluate_fitness()
                 initialized_model_IR12.update_individuals(1)
             return
 
         elif on_extinction == 'warn':
             with pytest.warns(RuntimeWarning, match="no offspring"):
-                initialized_model_IR12.evaluate_fitness(1)
+                initialized_model_IR12.evaluate_fitness()
                 initialized_model_IR12.update_individuals(1)
                 current = get_pop_subset()
-                initialized_model_IR12.evaluate_fitness(1)
+                initialized_model_IR12.evaluate_fitness()
                 initialized_model_IR12.update_individuals(1)
                 next = get_pop_subset()
 
         else:
-            initialized_model_IR12.evaluate_fitness(1)
+            initialized_model_IR12.evaluate_fitness()
             initialized_model_IR12.update_individuals(1)
             current = get_pop_subset()
-            initialized_model_IR12.evaluate_fitness(1)
+            initialized_model_IR12.evaluate_fitness()
             initialized_model_IR12.update_individuals(1)
             next = get_pop_subset()
 
@@ -354,7 +354,7 @@ def test_taxon_def(self, grid, trait_funcs, taxon_def, num_gen=10, dt=1):
 
         dfs = []
         for step in range(num_gen):
-            model.evaluate_fitness(dt)
+            model.evaluate_fitness()
             dfs.append(model.to_dataframe())
             model.update_individuals(dt)
 
@@ -375,5 +375,5 @@ def test_high_abundance_warning(self, grid, trait_funcs, num_gen=2, dt=1):
 
         with pytest.warns(RuntimeWarning, match="Large number of individuals generated"):
             for step in range(num_gen):
-                model.evaluate_fitness(dt)
+                model.evaluate_fitness()
                 model.update_individuals(dt)

adascape/tests/test_fastscape_ext.py

@@ -71,7 +71,7 @@ def test_speciation(speciation, specIR12_process):
     if speciation == 'IR12':
         spec = copy.deepcopy(specIR12_process)
     spec.initialize()
-    spec.run_step(1)
+    spec.run_step()
 
     assert spec.abundance == 10
     np.testing.assert_equal(spec._get_taxon_id(), np.ones(10))