-
Notifications
You must be signed in to change notification settings - Fork 37
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Merge branch 'main' into test_unitful_cellarea
- Loading branch information
Showing
31 changed files
with
377 additions
and
225 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -1,7 +1,7 @@ | ||
| name = "Rasters" | ||
| uuid = "a3a2b9e3-a471-40c9-b274-f788e487c689" | ||
| authors = ["Rafael Schouten <[email protected]>"] | ||
| version = "0.12.0" | ||
| version = "0.12.1" | ||
|
|
||
| [deps] | ||
| Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" | ||
|
|
@@ -25,6 +25,7 @@ ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca" | |
| RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" | ||
| Reexport = "189a3867-3050-52da-a836-e630ba90ab69" | ||
| Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46" | ||
| Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" | ||
|
|
||
| [weakdeps] | ||
| ArchGDAL = "c9ce4bd3-c3d5-55b8-8973-c0e20141b8c3" | ||
|
|
@@ -58,7 +59,7 @@ CommonDataModel = "0.2.3, 0.3" | |
| ConstructionBase = "1" | ||
| CoordinateTransformations = "0.6.2" | ||
| DataFrames = "1" | ||
| DimensionalData = "0.28.2" | ||
| DimensionalData = "0.29.4" | ||
| DiskArrays = "0.3, 0.4" | ||
| Extents = "0.1" | ||
| FillArrays = "0.12, 0.13, 1" | ||
|
|
||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -1,54 +1,119 @@ | ||
| Load occurrences for the Mountain Pygmy Possum using GBIF.jl | ||
| # Species distribution modelling workflow | ||
|
|
||
| ## Load GBIF | ||
| This example shows a full Species distribution modelling workflow, from loading data, to cleaning it, to fitting an ensemble and generating predictions. | ||
|
|
||
| ````@example gbif | ||
| using Rasters, GBIF2 | ||
| using RasterDataSources | ||
| const RS = Rasters | ||
| ```` | ||
| It uses GBIF and WorldClim data, which are common datasets in ecology. | ||
|
|
||
| ## Load Rasters, ArchGDAL, RasterDataSources and GBIF | ||
| The GBIF2 library is used to download occurrence data, RasterDataSources to conveniently access Bioclim data. ArchGDAL is necessary to load in the Bioclim data. | ||
|
|
||
| ````@example gbif | ||
| records = GBIF2.occurrence_search("Burramys parvus"; limit=300) | ||
| using Rasters, GBIF2 | ||
| using RasterDataSources, ArchGDAL | ||
| ```` | ||
|
|
||
| ## Extract coordinates | ||
|
|
||
| Extract the longitude/latitude value to a `Vector` of points | ||
| (a `Tuple` counts as a `(x, y)` point in GeoInterface.jl): | ||
| Load occurrences for the Mountain Pygmy Possum using GBIF.jl | ||
|
|
||
| ````@example gbif | ||
| coords = [(r.decimalLongitude, r.decimalLatitude) for r in records if !ismissing(r.decimalLatitude)] | ||
| records = GBIF2.occurrence_search("Burramys parvus"; limit=300) | ||
| ```` | ||
|
|
||
| ## Get layer / Band | ||
| Get BioClim layers and subset to south-east Australia | ||
| ## Get Bioclimatic variables | ||
| Get BioClim layers and subset to south-east Australia. | ||
| The first time this is run, this will automatically download and save the files. | ||
|
|
||
| ````@example gbif | ||
| A = RasterStack(WorldClim{BioClim}, (1, 3, 7, 12)) | ||
| se_aus = A[X(138 .. 155), Y(-40 .. -25), RS.Band(1)] | ||
| se_aus = A[X(138 .. 155), Y(-40 .. -25), Band(1)] | ||
| ```` | ||
| Plot BioClim predictors and scatter occurrence points on all subplots | ||
|
|
||
| ````@example gbif | ||
| using Plots | ||
| p = plot(se_aus); | ||
| kw = (legend=:none, opacity=0.5, markershape=:cross, markercolor=:black) | ||
| foreach(i -> scatter!(p, coords; subplot=i, kw...), 1:4) | ||
| # The coordinates from the gbif table | ||
| coords = collect(skipmissing(records.geometry)) | ||
| using CairoMakie | ||
| p = Rasters.rplot(se_aus); | ||
| for ax in p.content | ||
| if ax isa Axis | ||
| scatter!(ax, coords; alpha=0.5, marker='+', color=:black, markersize = 20) | ||
| end | ||
| end | ||
| p | ||
| ```` | ||
|
|
||
| Then extract predictor variables and write to CSV. | ||
| ## Extract bioclim variables at occurrence points | ||
| Then extract predictor variables and write to CSV. Use the skipmissing keyword to exclude both missing coordinates and coordinates with missing values in the RasterStack. | ||
|
|
||
| ````@example gbif | ||
| using CSV | ||
| predictors = collect(extract(se_aus, coords)) | ||
| CSV.write("burramys_parvus_predictors.csv", predictors) | ||
| presences = extract(se_aus, coords, skipmissing = true) | ||
| CSV.write("burramys_parvus_predictors.csv", presences) | ||
| ```` | ||
|
|
||
| Or convert them to a `DataFrame`: | ||
|
|
||
| ````@example gbif | ||
| using DataFrames | ||
| df = DataFrame(predictors) | ||
| df = DataFrame(presences) | ||
| df[1:5,:] | ||
| ```` | ||
| ```` | ||
|
|
||
| ## Sample background points | ||
| Next, sample random background points in the Raster. Rasters has a StatsBase extension to make this very straightforward. The syntax and output of `Rasters.sample` is very similar to that of `extract`. | ||
|
|
||
| ````@example gbif | ||
| using StatsBase | ||
| background = Rasters.sample(se_aus, 500, skipmissing = true) | ||
| ```` | ||
|
|
||
| ## Fit a statistical ensemble | ||
| In this example, we will [SpeciesDistributionModels.jl](https://github.com/tiemvanderdeure/SpeciesDistributionModels.jl) to fit a statistical ensemble to the occurrence and background data. | ||
|
|
||
| First we need to load the models. SDM.jl integrates with MLJ - see the [model browser](https://juliaai.github.io/MLJ.jl/dev/model_browser/#Classification) for what models are available. | ||
|
|
||
| ````@example gbif | ||
| import Maxnet: MaxnetBinaryClassifier | ||
| import MLJGLMInterface: LinearBinaryClassifier | ||
| # define the models in the ensemble | ||
| models = ( | ||
| maxnet = MaxnetBinaryClassifier(), | ||
| maxnet2 = MaxnetBinaryClassifier(features = "lq"), | ||
| glm = LinearBinaryClassifier() | ||
| ) | ||
| ```` | ||
|
|
||
| Next, format the data using `sdmdata`. To test how rigurous our models are, we will use 3-fold cross-validation. | ||
|
|
||
| ````@example gbif | ||
| using SpeciesDistributionModels | ||
| const SDM = SpeciesDistributionModels | ||
| data = sdmdata(presences, background; resampler = CV(; nfolds = 3)) | ||
| ```` | ||
|
|
||
| Now, fit the ensemble, passing the data object and the `NamedTuple` of models! | ||
|
|
||
| ````@example gbif | ||
| ensemble = sdm(data, models) | ||
| ```` | ||
|
|
||
| Use SDM.jl's evaluate function to see how this ensemble performs. | ||
|
|
||
| ````@example gbif | ||
| SDM.evaluate(ensemble) | ||
| ```` | ||
|
|
||
| Not too bad! | ||
|
|
||
| ## Make predictions of climatic suitability | ||
| Use the ensemble to | ||
|
|
||
| ````@example gbif | ||
| suitability = SDM.predict(ensemble, se_aus, reducer = mean) | ||
| ```` | ||
|
|
||
| And let's see what that looks like | ||
|
|
||
| ````@example gbif | ||
| plot(suitability, colorrange = (0,1)) | ||
| ```` |
Oops, something went wrong.