From cdb90358afc39a95cfa267bab2d02a34084705f1 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Thu, 30 Dec 2021 16:58:15 -0500 Subject: [PATCH 01/31] :art: temporal change prototype --- src/NeutralLandscapes.jl | 7 ++++- src/update.jl | 56 ++++++++++++++++++++++++++++++++++++++++ 2 files changed, 62 insertions(+), 1 deletion(-) create mode 100644 src/update.jl diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index bb8a7d6..6da5899 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -1,7 +1,7 @@ module NeutralLandscapes import NaNMath -using StatsBase: sample +using StatsBase: sample, ZScoreTransform, fit, transform using Random: rand! using Statistics: quantile, mean using Distributions: Normal @@ -40,4 +40,9 @@ include("algorithms/planargradient.jl") include("algorithms/rectangularcluster.jl") include("algorithms/wavesurface.jl") +include("update.jl") +export NeutralLandscapeUpdater, update!, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater + end # module + + diff --git a/src/update.jl b/src/update.jl new file mode 100644 index 0000000..ba208d2 --- /dev/null +++ b/src/update.jl @@ -0,0 +1,56 @@ +""" + This file contains methods for simulating change in a + landscape with expected amounts of spatial and temporal + autocorrelation. +""" + +""" + NeutralLandscapeUpdater is an abstract type for methods + for updating a landscape matrix +""" +abstract type NeutralLandscapeUpdater end + +""" + All `NeutralLandscapeUpdater`s have a field `rate` + which defines the expected change in any cell per timestep. +""" +rate(up::NeutralLandscapeUpdater) = up.rate + +generator(up::NeutralLandscapeUpdater) = up.generator + +function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} + _update(updater, mat) +end + + +@kwdef struct SpatiallyAutocorrelatedUpdater{G,R} <: NeutralLandscapeUpdater + generator::G = DiamondSquare(0.5) + rate::R = 0.1 +end +function _update(sau::SpatiallyAutocorrelatedUpdater, mat) + change = rand(generator(sau), size(mat)) + delta = rate(sau) .* transform(fit(ZScoreTransform, change), change) + mat .+ delta +end + + +@kwdef struct SpatiotemporallyAutocorrelatedUpdater{G,R,T} <: NeutralLandscapeUpdater + generator::G = DiamondSquare(0.1) + rate::R = 0.1 + temporalautocorrelation = 0.1 +end +function _update(sau::SpatiallyAutocorrelatedUpdater, mat) + +end + + + + +""" + example of how its used + +""" +#env = rand(MidpointDisplacement(0.9), (250,250)) +#sau = SpatiallyAutocorrelatedUpdater() +#new = update!(sau, env) +#plot(heatmap.([env, new])...) \ No newline at end of file From 3996277ae1c2e4bbacb557feb42ec546bc19b922 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Fri, 31 Dec 2021 14:40:12 -0500 Subject: [PATCH 02/31] spatiotemporal change --- src/NeutralLandscapes.jl | 28 ++++++++++--------- src/{algorithms => makers}/diamondsquare.jl | 0 src/{algorithms => makers}/discretevoronoi.jl | 0 .../distancegradient.jl | 0 src/{algorithms => makers}/edgegradient.jl | 0 src/{algorithms => makers}/nncluster.jl | 0 src/{algorithms => makers}/nnelement.jl | 0 src/{algorithms => makers}/nogradient.jl | 0 src/{algorithms => makers}/perlinnoise.jl | 0 src/{algorithms => makers}/planargradient.jl | 0 .../rectangularcluster.jl | 0 src/{algorithms => makers}/wavesurface.jl | 0 src/updaters/temporal.jl | 19 +++++++++++++ src/{ => updaters}/update.jl | 18 ++++++------ 14 files changed, 44 insertions(+), 21 deletions(-) rename src/{algorithms => makers}/diamondsquare.jl (100%) rename src/{algorithms => makers}/discretevoronoi.jl (100%) rename src/{algorithms => makers}/distancegradient.jl (100%) rename src/{algorithms => makers}/edgegradient.jl (100%) rename src/{algorithms => makers}/nncluster.jl (100%) rename src/{algorithms => makers}/nnelement.jl (100%) rename src/{algorithms => makers}/nogradient.jl (100%) rename src/{algorithms => makers}/perlinnoise.jl (100%) rename src/{algorithms => makers}/planargradient.jl (100%) rename src/{algorithms => makers}/rectangularcluster.jl (100%) rename src/{algorithms => makers}/wavesurface.jl (100%) create mode 100644 src/updaters/temporal.jl rename src/{ => updaters}/update.jl (71%) diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index 6da5899..434186a 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -28,19 +28,21 @@ export WaveSurface include("landscape.jl") include("classify.jl") -include("algorithms/diamondsquare.jl") -include("algorithms/discretevoronoi.jl") -include("algorithms/distancegradient.jl") -include("algorithms/edgegradient.jl") -include("algorithms/nncluster.jl") -include("algorithms/nnelement.jl") -include("algorithms/nogradient.jl") -include("algorithms/perlinnoise.jl") -include("algorithms/planargradient.jl") -include("algorithms/rectangularcluster.jl") -include("algorithms/wavesurface.jl") - -include("update.jl") +include("makers/diamondsquare.jl") +include("makers/discretevoronoi.jl") +include("makers/distancegradient.jl") +include("makers/edgegradient.jl") +include("makers/nncluster.jl") +include("makers/nnelement.jl") +include("makers/nogradient.jl") +include("makers/perlinnoise.jl") +include("makers/planargradient.jl") +include("makers/rectangularcluster.jl") +include("makers/wavesurface.jl") + +include("updaters/update.jl") +include("updaters/temporal.jl") +export TemporalUpdater, BrownianMotion, OhrsteinUhlenbeck export NeutralLandscapeUpdater, update!, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater end # module diff --git a/src/algorithms/diamondsquare.jl b/src/makers/diamondsquare.jl similarity index 100% rename from src/algorithms/diamondsquare.jl rename to src/makers/diamondsquare.jl diff --git a/src/algorithms/discretevoronoi.jl b/src/makers/discretevoronoi.jl similarity index 100% rename from src/algorithms/discretevoronoi.jl rename to src/makers/discretevoronoi.jl diff --git a/src/algorithms/distancegradient.jl b/src/makers/distancegradient.jl similarity index 100% rename from src/algorithms/distancegradient.jl rename to src/makers/distancegradient.jl diff --git a/src/algorithms/edgegradient.jl b/src/makers/edgegradient.jl similarity index 100% rename from src/algorithms/edgegradient.jl rename to src/makers/edgegradient.jl diff --git a/src/algorithms/nncluster.jl b/src/makers/nncluster.jl similarity index 100% rename from src/algorithms/nncluster.jl rename to src/makers/nncluster.jl diff --git a/src/algorithms/nnelement.jl b/src/makers/nnelement.jl similarity index 100% rename from src/algorithms/nnelement.jl rename to src/makers/nnelement.jl diff --git a/src/algorithms/nogradient.jl b/src/makers/nogradient.jl similarity index 100% rename from src/algorithms/nogradient.jl rename to src/makers/nogradient.jl diff --git a/src/algorithms/perlinnoise.jl b/src/makers/perlinnoise.jl similarity index 100% rename from src/algorithms/perlinnoise.jl rename to src/makers/perlinnoise.jl diff --git a/src/algorithms/planargradient.jl b/src/makers/planargradient.jl similarity index 100% rename from src/algorithms/planargradient.jl rename to src/makers/planargradient.jl diff --git a/src/algorithms/rectangularcluster.jl b/src/makers/rectangularcluster.jl similarity index 100% rename from src/algorithms/rectangularcluster.jl rename to src/makers/rectangularcluster.jl diff --git a/src/algorithms/wavesurface.jl b/src/makers/wavesurface.jl similarity index 100% rename from src/algorithms/wavesurface.jl rename to src/makers/wavesurface.jl diff --git a/src/updaters/temporal.jl b/src/updaters/temporal.jl new file mode 100644 index 0000000..d5811b9 --- /dev/null +++ b/src/updaters/temporal.jl @@ -0,0 +1,19 @@ +abstract type TemporalUpdater end + +update!(tu::T, x) where {T<:TemporalUpdater} = _getchange(tu, x) + +@kwdef struct BrownianMotion{S} <: TemporalUpdater + σ::S = 0.1 +end + +_getchange(bm::BrownianMotion, x) = rand(Normal(bm.σ)) + + +@kwdef struct OhrsteinUhlenbeck{S,T,F} <: TemporalUpdater + θ::T = 1. + σ::S = 0.1 + dt::F = 0.01 +end + +_getchange(ou::OhrsteinUhlenbeck, x) = ou.dt*(-ou.θ*x + rand(Normal(ou.σ))) + diff --git a/src/update.jl b/src/updaters/update.jl similarity index 71% rename from src/update.jl rename to src/updaters/update.jl index ba208d2..8fe83f5 100644 --- a/src/update.jl +++ b/src/updaters/update.jl @@ -21,26 +21,27 @@ generator(up::NeutralLandscapeUpdater) = up.generator function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) end - - @kwdef struct SpatiallyAutocorrelatedUpdater{G,R} <: NeutralLandscapeUpdater generator::G = DiamondSquare(0.5) rate::R = 0.1 end function _update(sau::SpatiallyAutocorrelatedUpdater, mat) change = rand(generator(sau), size(mat)) - delta = rate(sau) .* transform(fit(ZScoreTransform, change), change) + delta = rate(sau) .+ rate(sau) .* transform(fit(ZScoreTransform, change), change) mat .+ delta end @kwdef struct SpatiotemporallyAutocorrelatedUpdater{G,R,T} <: NeutralLandscapeUpdater generator::G = DiamondSquare(0.1) - rate::R = 0.1 - temporalautocorrelation = 0.1 + rate::R = 0.1 + temporalupdater::T = BrownianMotion() end -function _update(sau::SpatiallyAutocorrelatedUpdater, mat) - +function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) + change = rand(generator(stau), size(mat)) + temporalshift = broadcast(x->update!(stau.temporalupdater, x), mat) + delta = rate(stau) .+ temporalshift .* transform(fit(ZScoreTransform, change), change) + mat .+ delta end @@ -53,4 +54,5 @@ end #env = rand(MidpointDisplacement(0.9), (250,250)) #sau = SpatiallyAutocorrelatedUpdater() #new = update!(sau, env) -#plot(heatmap.([env, new])...) \ No newline at end of file +#plot(heatmap.([env, new])...) + From e04e9db5f610dc419a2b75df41d46f84dc893a2b Mon Sep 17 00:00:00 2001 From: michael catchen Date: Fri, 31 Dec 2021 14:41:24 -0500 Subject: [PATCH 03/31] update returns, not inplace --- src/NeutralLandscapes.jl | 2 +- src/updaters/temporal.jl | 2 +- src/updaters/update.jl | 4 ++-- 3 files changed, 4 insertions(+), 4 deletions(-) diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index 434186a..e25dc8a 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -43,7 +43,7 @@ include("makers/wavesurface.jl") include("updaters/update.jl") include("updaters/temporal.jl") export TemporalUpdater, BrownianMotion, OhrsteinUhlenbeck -export NeutralLandscapeUpdater, update!, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater +export NeutralLandscapeUpdater, update, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater end # module diff --git a/src/updaters/temporal.jl b/src/updaters/temporal.jl index d5811b9..cbf296c 100644 --- a/src/updaters/temporal.jl +++ b/src/updaters/temporal.jl @@ -1,6 +1,6 @@ abstract type TemporalUpdater end -update!(tu::T, x) where {T<:TemporalUpdater} = _getchange(tu, x) +update(tu::T, x) where {T<:TemporalUpdater} = _getchange(tu, x) @kwdef struct BrownianMotion{S} <: TemporalUpdater σ::S = 0.1 diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 8fe83f5..55a2660 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -18,7 +18,7 @@ rate(up::NeutralLandscapeUpdater) = up.rate generator(up::NeutralLandscapeUpdater) = up.generator -function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} +function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) end @kwdef struct SpatiallyAutocorrelatedUpdater{G,R} <: NeutralLandscapeUpdater @@ -39,7 +39,7 @@ end end function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) change = rand(generator(stau), size(mat)) - temporalshift = broadcast(x->update!(stau.temporalupdater, x), mat) + temporalshift = broadcast(x->update(stau.temporalupdater, x), mat) delta = rate(stau) .+ temporalshift .* transform(fit(ZScoreTransform, change), change) mat .+ delta end From 9f40c87c2df0daced9ecf4c1dbb3dfd30680b9e0 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Fri, 31 Dec 2021 14:48:05 -0500 Subject: [PATCH 04/31] :bug: spatiallytemporal update fix --- src/updaters/update.jl | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 55a2660..3b4a51d 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -40,7 +40,8 @@ end function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) change = rand(generator(stau), size(mat)) temporalshift = broadcast(x->update(stau.temporalupdater, x), mat) - delta = rate(stau) .+ temporalshift .* transform(fit(ZScoreTransform, change), change) + z = transform(fit(ZScoreTransform, change), change) + delta = rate(stau) .+ (rate(stau) .* z) .+ (temporalshift .* z) mat .+ delta end From 78b3907f1eae5be2062e224909c7d7c0de61e436 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Fri, 31 Dec 2021 15:59:13 -0500 Subject: [PATCH 05/31] :wrench: folder change to a new pr --- src/NeutralLandscapes.jl | 22 ++-- src/makers/diamondsquare.jl | 206 ------------------------------- src/makers/discretevoronoi.jl | 21 ---- src/makers/distancegradient.jl | 24 ---- src/makers/edgegradient.jl | 21 ---- src/makers/nncluster.jl | 33 ----- src/makers/nnelement.jl | 43 ------- src/makers/nogradient.jl | 11 -- src/makers/perlinnoise.jl | 136 -------------------- src/makers/planargradient.jl | 24 ---- src/makers/rectangularcluster.jl | 29 ----- src/makers/wavesurface.jl | 22 ---- 12 files changed, 11 insertions(+), 581 deletions(-) delete mode 100644 src/makers/diamondsquare.jl delete mode 100644 src/makers/discretevoronoi.jl delete mode 100644 src/makers/distancegradient.jl delete mode 100644 src/makers/edgegradient.jl delete mode 100644 src/makers/nncluster.jl delete mode 100644 src/makers/nnelement.jl delete mode 100644 src/makers/nogradient.jl delete mode 100644 src/makers/perlinnoise.jl delete mode 100644 src/makers/planargradient.jl delete mode 100644 src/makers/rectangularcluster.jl delete mode 100644 src/makers/wavesurface.jl diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index e25dc8a..a3c6274 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -28,17 +28,17 @@ export WaveSurface include("landscape.jl") include("classify.jl") -include("makers/diamondsquare.jl") -include("makers/discretevoronoi.jl") -include("makers/distancegradient.jl") -include("makers/edgegradient.jl") -include("makers/nncluster.jl") -include("makers/nnelement.jl") -include("makers/nogradient.jl") -include("makers/perlinnoise.jl") -include("makers/planargradient.jl") -include("makers/rectangularcluster.jl") -include("makers/wavesurface.jl") +include("algorithms/diamondsquare.jl") +include("algorithms/discretevoronoi.jl") +include("algorithms/distancegradient.jl") +include("algorithms/edgegradient.jl") +include("algorithms/nncluster.jl") +include("algorithms/nnelement.jl") +include("algorithms/nogradient.jl") +include("algorithms/perlinnoise.jl") +include("algorithms/planargradient.jl") +include("algorithms/rectangularcluster.jl") +include("algorithms/wavesurface.jl") include("updaters/update.jl") include("updaters/temporal.jl") diff --git a/src/makers/diamondsquare.jl b/src/makers/diamondsquare.jl deleted file mode 100644 index a35506e..0000000 --- a/src/makers/diamondsquare.jl +++ /dev/null @@ -1,206 +0,0 @@ -""" - DiamondSquare <: NeutralLandscapeMaker - - DiamondSquare(; H = 0.5) - DiamondSquare(H) - -This type generates a neutral landscape using the diamond-squares -algorithm, which produces fractals with variable spatial autocorrelation. - -https://en.wikipedia.org/wiki/Diamond-square_algorithm - -The algorithm is named diamond-square because it is an iterative procedure of -"diamond" and "square" steps. - - -The degree of spatial autocorrelation is controlled by a parameter `H`, -which varies from 0.0 (low autocorrelation) to 1.0 (high autocorrelation) --- -note this is non-inclusive and H = 0 and H = 1 will not behave as expected. -The result of the diamond-square algorithm is a fractal with dimension D = 2 + H. - -A similar algorithm, midpoint-displacement, functions almost -identically, except that in DiamondSquare, the square step interpolates -edge midpoints from the nearest two corners and the square's center, where as -midpoint-displacement only interpolates from the nearest corners (see `MidpointDisplacement`). - -""" -@kwdef struct DiamondSquare <: NeutralLandscapeMaker - H::Float64 = 0.5 - function DiamondSquare(H::T) where {T <: Real} - @assert 0 <= H < 1 - new(H) - end -end - -""" - MidpointDisplacement <: NeutralLandscapeMaker - - MidpointDisplacement(; H = 0.5) - -Creates a midpoint-displacement algorithm object `MidpointDisplacement`. -The degree of spatial autocorrelation is controlled by a parameter `H`, -which varies from 0.0 (low autocorrelation) to 1.0 (high autocorrelation) --- -note this is non-inclusive and H = 0 and H = 1 will not behave as expected. - -A similar algorithm, diamond-square, functions almost -identically, except that in diamond-square, the square step interpolates -edge midpoints from the nearest two corners and the square's center, where as -`MidpointDisplacement` only interpolates from the nearest corners (see `DiamondSquare`). -""" -@kwdef struct MidpointDisplacement <: NeutralLandscapeMaker - H::Float64 = 0.5 - function MidpointDisplacement(H::T) where {T <: Real} - @assert 0 <= H < 1 - new(H) - end -end - -# Check if `mat` is the right size. If mat is not the correct size (DiamondSquare -# can only run on a lattice of size NxN where N = (2^n)+1 for integer n), -# allocates the smallest lattice large enough to contain `mat` that can run -# DiamondSquare. -function _landscape!(mat, alg::Union{DiamondSquare, MidpointDisplacement}; kw...) where {IT <: Integer} - - rightSize::Bool = _isPowerOfTwo(size(mat)[1]-1) && _isPowerOfTwo(size(mat)[2]-1) - latticeSize::Int = size(mat)[1] - - dsMat = mat - if !rightSize - dim1, dim2 = size(mat) - smallestContainingLattice::Int = 2^ceil(log2(max(dim1, dim2))) + 1 - dsMat = zeros(smallestContainingLattice, smallestContainingLattice) - end - _diamondsquare!(dsMat, alg) - - mat .= dsMat[1:size(mat)[1], 1:size(mat)[2]] -end - -# Runs the diamond-square algorithm on a matrix `mat` of size -# `NxN`, where `N=(2^n)+1` for some integer `n`, i.e (N=5,9,17,33,65) - -# Diamond-square is an iterative procedure, where the lattice is divided -# into subsquares in subsequent rounds. At each round, the subsquares shrink in size, -# as previously uninitialized values in the lattice are interpolated as a mean of nearby points plus random displacement. -# As the rounds increase, the magnitude of this displacement decreases. This creates spatioautocorrelation, which is controlled -# by a single parameter `H` which varies between `0` (no autocorrelation) and `1` (high autocorrelation) -function _diamondsquare!(mat, alg) - latticeSize = size(mat)[1] - numberOfRounds::Int = log2(latticeSize-1) - _initializeDiamondSquare!(mat, alg) - - for round in 0:(numberOfRounds-1) # counting from 0 saves us a headache later - subsquareSideLength::Int = 2^(numberOfRounds-(round)) - numberOfSubsquaresPerAxis::Int = ((latticeSize-1) / subsquareSideLength)-1 - - for x in 0:numberOfSubsquaresPerAxis # iterate over the subsquares within the lattice at this side length - for y in 0:numberOfSubsquaresPerAxis - subsquareCorners = _subsquareCornerCoordinates(x,y,subsquareSideLength) - - _diamond!(mat, alg, round, subsquareCorners) - _square!(mat, alg, round, subsquareCorners) - end - end - end -end - -# Initialize's the `DiamondSquare` algorithm by displacing the four corners of the -# lattice using `displace`, scaled by the algorithm's autocorrelation `H`. -function _initializeDiamondSquare!(mat, alg) - latticeSize = size(mat)[1] - corners = _subsquareCornerCoordinates(0,0, latticeSize-1) - for mp in corners - mat[mp...] = _displace(alg.H, 1) - end -end - -# Returns the coordinates for the corners of the subsquare (x,y) given a side-length `sideLength`. -function _subsquareCornerCoordinates(x::Int, y::Int, sideLength::Int) - corners = [1 .+ sideLength.*i for i in [(x,y), (x+1, y), (x, y+1), (x+1, y+1)]] -end - -# Runs the diamond step of the `DiamondSquare` algorithm on the square defined by -# `corners` on the matrix `mat`. The center of the square is interpolated from the -# four corners, and is displaced. The displacement is drawn according to `alg.H` and round using `displace` -function _diamond!(mat, alg, round::Int, corners::AbstractVector{Tuple{Int, Int}}) - centerPt = _centerCoordinate(corners) - mat[centerPt...] = _interpolate(mat, corners) + _displace(alg.H, round) -end - -# Runs the square step of the `DiamondSquare` algorithm on the square defined -# by `corners` on the matrix `mat`. The midpoint of each edge of this square is interpolated -# by computing the mean value of the two corners on the edge and the center of the square, and the -# displacing it. The displacement is drawn according to `alg.H` and round using `displace` -function _square!(mat, alg::DiamondSquare, round::Int, corners::AbstractVector{Tuple{Int, Int}}) - bottomLeft,bottomRight,topLeft,topRight = corners - leftEdge, bottomEdge, topEdge, rightEdge = _edgeMidpointCoordinates(corners) - centerPoint = _centerCoordinate(corners) - - mat[leftEdge...] = _interpolate(mat, [topLeft,bottomLeft,centerPoint]) + _displace(alg.H, round) - mat[bottomEdge...] = _interpolate(mat, [bottomLeft,bottomRight,centerPoint]) + _displace(alg.H, round) - mat[topEdge...] = _interpolate(mat, [topLeft,topRight,centerPoint]) + _displace(alg.H, round) - mat[rightEdge...] = _interpolate(mat, [topRight,bottomRight,centerPoint]) + _displace(alg.H, round) -end - -# Runs the square step of the `MidpointDisplacement` algorithm on the square defined -# by `corners` on the matrix `mat`. The midpoint of each edge of this square is interpolated -# by computing the mean value of the two corners on the edge and the center of the square, and the -# displacing it. The displacement is drawn according to `alg.H` and round using `displace` -function _square!(mat, alg::MidpointDisplacement, round::Int, corners::AbstractVector{Tuple{Int, Int}}) - bottomLeft,bottomRight,topLeft,topRight = corners - leftEdge, bottomEdge, topEdge, rightEdge = _edgeMidpointCoordinates(corners) - mat[leftEdge...] = _interpolate(mat, [topLeft,bottomLeft]) + _displace(alg.H, round) - mat[bottomEdge...] = _interpolate(mat, [bottomLeft,bottomRight]) + _displace(alg.H, round) - mat[topEdge...] = _interpolate(mat, [topLeft,topRight]) + _displace(alg.H, round) - mat[rightEdge...] = _interpolate(mat, [topRight,bottomRight]) + _displace(alg.H, round) -end - -# Computes the mean of a set of points, represented as a list of indicies to a matrix `mat`. -function _interpolate(mat, points::AbstractVector{Tuple{Int,Int}}) - return mean(mat[pt...] for pt in points) -end - -# `displace` produces a random value as a function of `H`, which is the -# autocorrelation parameter used in `DiamondSquare` and must be between `0` -# and `1`, and `round` which describes the current tiling size for the -# DiamondSquare() algorithm. - -# Random value are drawn from a Gaussian distribution using `Distribution.Normal` -# The standard deviation of this Gaussian, σ, is set to (1/2)^(round*H), which will -# move from 1.0 to 0 as `round` increases. -function _displace(H::Float64, round::Int) - σ = 0.5^(round*H) - return(rand(Normal(0, σ))) -end - -# Returns the center coordinate for a square defined by `corners` for the -# `DiamondSquare` algorithm. -function _centerCoordinate(corners::AbstractVector{Tuple{Int,Int}}) - bottomLeft,bottomRight,topLeft,topRight = corners - centerX::Int = (_xcoord(bottomLeft)+_xcoord(bottomRight)) ÷ 2 - centerY::Int = (_ycoord(topRight)+_ycoord(bottomRight)) ÷ 2 - return (centerX, centerY) -end - -# Returns the x-coordinate from a lattice coordinate `pt`. -_xcoord(pt::Tuple{Int,Int}) = pt[1] -# Returns the y-coordinate from a lattice coordinate `pt`. -_ycoord(pt::Tuple{Int,Int}) = pt[2] - -# Returns an array of midpoints for a square defined by `corners` for the `DiamondSquare` algorithm. -function _edgeMidpointCoordinates(corners::AbstractVector{Tuple{Int,Int}}) - # bottom left, bottom right, top left, top right - bottomLeft,bottomRight,topLeft,topRight = corners - - leftEdgeMidpoint::Tuple{Int,Int} = (_xcoord(bottomLeft), (_ycoord(bottomLeft)+_ycoord(topLeft))÷2 ) - bottomEdgeMidpoint::Tuple{Int,Int} = ( (_xcoord(bottomLeft)+ _xcoord(bottomRight))÷2, _ycoord(bottomLeft) ) - topEdgeMidpoint::Tuple{Int,Int} = ( (_xcoord(topLeft)+_xcoord(topRight))÷2, _ycoord(topLeft)) - rightEdgeMidpoint::Tuple{Int,Int} = ( _xcoord(bottomRight), (_ycoord(bottomRight)+_ycoord(topRight))÷2) - - edgeMidpoints = [leftEdgeMidpoint, bottomEdgeMidpoint, topEdgeMidpoint, rightEdgeMidpoint] - return edgeMidpoints -end - -# Determines if `x`, an integer, can be expressed as `2^n`, where `n` is also an integer. -function _isPowerOfTwo(x::IT) where {IT <: Integer} - return x & (x-1) == 0 -end diff --git a/src/makers/discretevoronoi.jl b/src/makers/discretevoronoi.jl deleted file mode 100644 index 8155083..0000000 --- a/src/makers/discretevoronoi.jl +++ /dev/null @@ -1,21 +0,0 @@ -""" - DiscreteVoronoi <: NeutralLandscapeMaker - - DiscreteVoronoi(; n=3) - DiscreteVoronoi(n) - -This type provides a rasterization of a Voronoi-like diagram. -Assigns a value to each patch using a 1-NN algorithmm with `n` initial clusters. -It is a `NearestNeighborElement` algorithmm with `k` neighbors set to 1. -The default is to use three clusters. -""" -@kwdef struct DiscreteVoronoi <: NeutralLandscapeMaker - n::Int64 = 3 - function DiscreteVoronoi(n::Int64) - @assert n > 0 - new(n) - end -end - -_landscape!(mat, alg::DiscreteVoronoi) = _landscape!(mat, NearestNeighborElement(alg.n, 1)) - diff --git a/src/makers/distancegradient.jl b/src/makers/distancegradient.jl deleted file mode 100644 index d2e2b1d..0000000 --- a/src/makers/distancegradient.jl +++ /dev/null @@ -1,24 +0,0 @@ -""" - DistanceGradient <: NeutralLandscapeMaker - - DistanceGradient(; sources=[1]) - DistanceGradient(sources) - -The `sources` field is a `Vector{Integer}` of *linear* indices of the matrix, -from which the distance must be calculated. -""" -@kwdef struct DistanceGradient <: NeutralLandscapeMaker - sources::Vector{Integer} = [1] -end - -function _landscape!(mat, alg::DistanceGradient) - @assert maximum(alg.sources) <= length(mat) - @assert minimum(alg.sources) > 0 - coordinates = Matrix{Float64}(undef, (2, prod(size(mat)))) - for (i, p) in enumerate(Iterators.product(axes(mat)...)) - coordinates[1:2, i] .= p - end - tree = KDTree(coordinates[:,alg.sources]) - mat[:] .= nn(tree, coordinates)[2] - return mat -end diff --git a/src/makers/edgegradient.jl b/src/makers/edgegradient.jl deleted file mode 100644 index e0839f5..0000000 --- a/src/makers/edgegradient.jl +++ /dev/null @@ -1,21 +0,0 @@ -""" - EdgeGradient <: NeutralLandscapeMaker - - EdgeGradient(; direction=360rand()) - EdgeGradient(direction) - -This type is used to generate an edge gradient landscape, where values change -as a bilinear function of the *x* and *y* coordinates. The direction is -expressed as a floating point value, which will be in *[0,360]*. The inner -constructor takes the mod of the value passed and 360, so that a value that is -out of the correct interval will be corrected. -""" -@kwdef struct EdgeGradient <: NeutralLandscapeMaker - direction::Float64 = 360rand() - EdgeGradient(x::T) where {T <: Real} = new(mod(x, 360.0)) -end - -function _landscape!(mat, alg::EdgeGradient) where {IT <: Integer} - _landscape!(mat, PlanarGradient(alg.direction)) - mat .= -2.0abs.(0.5 .- mat) .+ 1.0 -end diff --git a/src/makers/nncluster.jl b/src/makers/nncluster.jl deleted file mode 100644 index a52af29..0000000 --- a/src/makers/nncluster.jl +++ /dev/null @@ -1,33 +0,0 @@ -""" - NearestNeighborCluster <: NeutralLandscapeMaker - - NearestNeighborCluster(; p=0.5, n=:rook) - NearestNeighborCluster(p, [n=:rook]) - -Create a random cluster nearest-neighbour neutral landscape model with -values ranging 0-1. `p` sets the density of original clusters, and `n` -sets the neighborhood for clustering (see `?label` for neighborhood options) -""" -@kwdef struct NearestNeighborCluster <: NeutralLandscapeMaker - p::Float64 = 0.5 - n::Symbol = :rook - function NearestNeighborCluster(p::Float64, n::Symbol = :rook) - @assert p > 0 - @assert n ∈ (:rook, :queen, :diagonal) - new(p,n) - end -end - -function _landscape!(mat, alg::NearestNeighborCluster) - _landscape!(mat, NoGradient()) - classify!(mat, [alg.p, 1 - alg.p]) - replace!(mat, 2.0 => NaN) - clusters, nClusters = label(mat, alg.n) - coordinates = _coordinatematrix(clusters) - sources = findall(!isnan, vec(clusters)) - tree = KDTree(coordinates[:,sources]) - clusters[:] .= clusters[sources[nn(tree, coordinates)[1]]] - randvals = rand(nClusters) - mat .= randvals[Int.(clusters)] - mat -end diff --git a/src/makers/nnelement.jl b/src/makers/nnelement.jl deleted file mode 100644 index a5ba8cf..0000000 --- a/src/makers/nnelement.jl +++ /dev/null @@ -1,43 +0,0 @@ -""" - NearestNeighborElement <: NeutralLandscapeMaker - - NearestNeighborElement(; n=3, k=1) - NearestNeighborElement(n, [k=1]) - -Assigns a value to each patch using a k-NN algorithmm with `n` initial clusters -and `k` neighbors. The default is to use three cluster and a single neighbor. -""" -@kwdef struct NearestNeighborElement <: NeutralLandscapeMaker - n::Int64 = 3 - k::Int64 = 1 - function NearestNeighborElement(n::Int64, k::Int64 = 1) - @assert n > 0 - @assert k > 0 - @assert k <= n - new(n,k) - end -end - -function _coordinatematrix(mat) - coordinates = Matrix{Float64}(undef, (2, prod(size(mat)))) - for (i, p) in enumerate(Iterators.product(axes(mat)...)) - coordinates[1:2, i] .= p - end - coordinates -end - -function _landscape!(mat, alg::NearestNeighborElement) - fill!(mat, 0.0) - clusters = sample(eachindex(mat), alg.n; replace=false) - for (i,n) in enumerate(clusters) - mat[n] = i - end - coordinates = _coordinatematrix(mat) - tree = KDTree(coordinates[:,clusters]) - if alg.k == 1 - mat[:] .= nn(tree, coordinates)[1] - else - mat[:] .= map(mean, knn(tree, coordinates, alg.k)[1]) - end - return mat -end diff --git a/src/makers/nogradient.jl b/src/makers/nogradient.jl deleted file mode 100644 index c05f396..0000000 --- a/src/makers/nogradient.jl +++ /dev/null @@ -1,11 +0,0 @@ -""" - NoGradient <: NeutralLandscapeMaker - - NoGradient() - -This type is used to generate a random landscape with no gradients -""" -struct NoGradient <: NeutralLandscapeMaker -end - -_landscape!(mat, ::NoGradient) = rand!(mat) diff --git a/src/makers/perlinnoise.jl b/src/makers/perlinnoise.jl deleted file mode 100644 index 3265682..0000000 --- a/src/makers/perlinnoise.jl +++ /dev/null @@ -1,136 +0,0 @@ -""" - PerlinNoise <: NeutralLandscapeMaker - - PerlinNoise(; kw...) - PerlinNoise(periods, [octaves=1, lacunarity=2, persistance=0.5, valley=:u]) - -Create a Perlin noise neutral landscape model with values ranging 0-1. - -# Keywords - -- `periods::Tuple{Int,Int}=(1,1)`: the number of periods of Perlin noise across row and - column dimensions for the first octave. -- `octaves::Int=1`: the number of octaves that will form the Perlin noise. -- `lacunarity::Int=2` : the rate at which the frequency of periods increases for each - octive. -- `persistance::Float64=0.5` : the rate at which the amplitude of periods decreases for each - octive. -- `valley::Symbol=`:u`: the kind of valley bottom that will be mimicked: `:u` produces - u-shaped valleys, `:v` produces v-shaped valleys, and `:-` produces flat bottomed - valleys. - -Note: This is a memory-intensive algorithm with some settings. Be careful using larger -prime numbers for `period` when also using a large array size, high lacuarity and/or many -octaves. Memory use scales with the lowest common multiple of `periods`. -""" -@kwdef struct PerlinNoise <: NeutralLandscapeMaker - periods::Tuple{Int,Int} = (1, 1) - octaves::Int = 1 - lacunarity::Int = 2 - persistance::Float64 = 0.5 - valley::Symbol = :u - function PerlinNoise(periods, octaves=1, lacunarity=2, persistance=0.5, valley=:u) - new(periods, octaves, lacunarity, persistance, valley) - end -end - -function _landscape!(A, alg::PerlinNoise) - # nRow must equal nCol so determine the dimension of the smallest square - dim = max(size(A)...) - # Check the dim is a multiple of each octives maximum number of periods and - # expand dim if needed - rperiodsmax = alg.periods[1] * alg.lacunarity^(alg.octaves - 1) - cperiodsmax = alg.periods[2] * alg.lacunarity^(alg.octaves - 1) - periodsmultiple = lcm(rperiodsmax, cperiodsmax) # lowest common multiple - if dim % periodsmultiple != 0 - dim = ceil(Int, dim / periodsmultiple) * periodsmultiple - end - - # Generate the Perlin noise - noise = zeros(eltype(A), (dim, dim)) - meshbuf = Array{eltype(A),3}(undef, dim, dim, 2) - nbufs = ntuple(_->Array{eltype(A),2}(undef, dim, dim), 4) - for octave in 0:alg.octaves-1 - octave_noise!(noise, meshbuf, nbufs, alg, octave, dim, dim) - end - - # Randomly extract the desired array size - noiseview = _view_from_square(noise, size(A)...) - - # Rescale the Perlin noise to mimic different kinds of valley bottoms - return if alg.valley == :u - A .= noiseview - elseif alg.valley == :v - A .= abs.(noiseview) - elseif alg.valley == :- - A .= noiseview.^2 - else - error("$(alg.valley) not recognised for `valley` use `:u`, `:v` or `:-`") - end -end - -function octave_noise!( - noise, mesh, (n11, n21, n12, n22), alg::PerlinNoise, octave, nrow, ncol -) - f(t) = @fastmath 6 * t ^ 5 - 15 * t ^ 4 + 10 * t ^ 3 # Wut - - # Mesh - rp, cp = alg.periods .* alg.lacunarity^(octave) - delta = (rp / nrow, cp / ncol) - ranges = range(0, rp-delta[1], length=nrow), range(0, cp-delta[2], length=ncol) - _mesh!(mesh, ranges...) - @fastmath mesh .%= 1 - - # Gradients - # This allocates, but the gradients size changes with octave so it needs to - # be a very smart in-place `repeat!` or some kind of generator, and the improvement - # may not be than large (~20%). - angles = 2pi .* rand(rp + 1, cp + 1) - @fastmath gradients = cat(cos.(angles), sin.(angles); dims=3) - d = (nrow ÷ rp, ncol ÷ cp) - grad = repeat(gradients, inner=[d[1], d[2], 1]) - - g111 = @view grad[1:ncol, 1:ncol, 1] - g211 = @view grad[end-nrow+1:end, 1:ncol, 1] - g121 = @view grad[1:ncol, end-ncol+1:end, 1] - g221 = @view grad[end-nrow+1:end, end-ncol+1:end, 1] - g112 = @view grad[1:ncol, 1:ncol, 2] - g212 = @view grad[end-nrow+1:end, 1:ncol, 2] - g122 = @view grad[1:ncol, end-ncol+1:end, 2] - g222 = @view grad[end-nrow+1:end, end-ncol+1:end, 2] - - # Ramps - m1 = @view mesh[:, :, 1] - m2 = @view mesh[:, :, 2] - n11 .= ((m1 .+ m2 ) .* g111 .+ (m1 .+ m2 ) .* g112) - n21 .= ((m1 .-1 .+ m2 ) .* g211 .+ (m1 .-1 .+ m2 ) .* g212) - n12 .= ((m1 .+ m2 .- 1) .* g121 .+ (m1 .+ m2 .- 1) .* g122) - n22 .= ((m1 .-1 .+ m2 .- 1) .* g221 .+ (m1 .-1 .+ m2 .- 1) .* g222) - - # Interpolation - mesh .= f.(mesh) - t1 = @view mesh[:, :, 1] - t2 = @view mesh[:, :, 2] - noise .+= sqrt(2) .* (alg.persistance ^ octave) .* - ((1 .- t2) .* (n11 .* (1 .- t1) .+ t1 .* n21) .+ - t2 .* (n12 .* (1 .- t1) .+ t1 .* n22)) - return noise -end - -function _mesh!(A, x, y) - for (i, ival) in enumerate(x), j in 1:length(y) - A[i, j, 1] = ival - end - for i in 1:length(x), (j, jval) in enumerate(y) - A[i, j, 2] = jval - end - return A -end - -function _view_from_square(source, nrow, ncol) - # Extract a portion of the array to match the dimensions - dim = size(source, 1) - startrow = rand(1:(dim - nrow + 1)) - startcol = rand(1:(dim - ncol + 1)) - return @view source[startrow:startrow + nrow - 1, startcol:startcol + ncol - 1] -end diff --git a/src/makers/planargradient.jl b/src/makers/planargradient.jl deleted file mode 100644 index 601b06b..0000000 --- a/src/makers/planargradient.jl +++ /dev/null @@ -1,24 +0,0 @@ -""" - PlanarGradient <: NeutralLandscapeMaker - - PlanarGradient(; direction=360rand()) - PlanarGradient(direction) - -This type is used to generate a planar gradient landscape, where values change -as a bilinear function of the *x* and *y* coordinates. The direction is -expressed as a floating point value, which will be in *[0,360]*. The inner -constructor takes the mod of the value passed and 360, so that a value that is -out of the correct interval will be corrected. -""" -@kwdef struct PlanarGradient <: NeutralLandscapeMaker - direction::Float64 = 360rand() - PlanarGradient(x::T) where {T <: Real} = new(mod(x, 360.0)) -end - -function _landscape!(mat, alg::PlanarGradient) where {IT <: Integer} - eastness = sin(deg2rad(alg.direction)) - southness = -1cos(deg2rad(alg.direction)) - rows, cols = axes(mat) - mat .= southness .* rows .+ (eastness .* cols)' - return mat -end diff --git a/src/makers/rectangularcluster.jl b/src/makers/rectangularcluster.jl deleted file mode 100644 index f2b6ff6..0000000 --- a/src/makers/rectangularcluster.jl +++ /dev/null @@ -1,29 +0,0 @@ -""" - RectangularCluster <: NeutralLandscapeMaker - - RectangularCluster(; minimum=2, maximum=4) - RectangularCluster(minimum, [maximum=4]) - -Fills the landscape with rectangles containing a random value. The size of each -rectangle/patch is between `minimum` and `maximum` (the two can be equal for a -fixed size rectangle). -""" -@kwdef struct RectangularCluster <: NeutralLandscapeMaker - minimum::Integer = 2 - maximum::Integer = 4 - function RectangularCluster(minimum::T, maximum::T=4) where {T <: Integer} - @assert 0 < minimum <= maximum - new(minimum, maximum) - end -end - -function _landscape!(mat, alg::RectangularCluster) - mat .= -1.0 - while minimum(mat) == -1.0 - width, height = rand(alg.minimum:alg.maximum, 2) - row = rand(1:(size(mat,1)-(width-1))) - col = rand(1:(size(mat,2)-(height-1))) - mat[row:(row+(width-1)) , col:(col+(height-1))] .= rand() - end - return mat -end diff --git a/src/makers/wavesurface.jl b/src/makers/wavesurface.jl deleted file mode 100644 index fa86d82..0000000 --- a/src/makers/wavesurface.jl +++ /dev/null @@ -1,22 +0,0 @@ -""" - WaveSurface <: NeutralLandscapeMaker - - WaveSurface(; direction=360rand(), periods=1) - WaveSurface(direction, [periods=1]) - -Creates a sinusoidal landscape with a `direction` and a number of `periods`. If -neither are specified, there will be a single period of random direction. -""" -@kwdef struct WaveSurface <: NeutralLandscapeMaker - direction::Float64 = 360rand() - periods::Int64 = 1 - function WaveSurface(direction::T, periods::K = 1) where {T <: Real, K <: Integer} - @assert periods >= 1 - new(mod(direction, 360.0), periods) - end -end - -function _landscape!(mat, alg::WaveSurface) - rand!(mat, PlanarGradient(alg.direction)) - mat .= sin.(mat .* (2π * alg.periods)) -end From 176f51eeba0fe6032942a797edf5783bc049b1c9 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Fri, 31 Dec 2021 16:00:27 -0500 Subject: [PATCH 06/31] does this revert work? --- src/algorithms/diamondsquare.jl | 206 +++++++++++++++++++++++++++ src/algorithms/discretevoronoi.jl | 21 +++ src/algorithms/distancegradient.jl | 24 ++++ src/algorithms/edgegradient.jl | 21 +++ src/algorithms/nncluster.jl | 33 +++++ src/algorithms/nnelement.jl | 43 ++++++ src/algorithms/nogradient.jl | 11 ++ src/algorithms/perlinnoise.jl | 136 ++++++++++++++++++ src/algorithms/planargradient.jl | 24 ++++ src/algorithms/rectangularcluster.jl | 29 ++++ src/algorithms/wavesurface.jl | 22 +++ 11 files changed, 570 insertions(+) create mode 100644 src/algorithms/diamondsquare.jl create mode 100644 src/algorithms/discretevoronoi.jl create mode 100644 src/algorithms/distancegradient.jl create mode 100644 src/algorithms/edgegradient.jl create mode 100644 src/algorithms/nncluster.jl create mode 100644 src/algorithms/nnelement.jl create mode 100644 src/algorithms/nogradient.jl create mode 100644 src/algorithms/perlinnoise.jl create mode 100644 src/algorithms/planargradient.jl create mode 100644 src/algorithms/rectangularcluster.jl create mode 100644 src/algorithms/wavesurface.jl diff --git a/src/algorithms/diamondsquare.jl b/src/algorithms/diamondsquare.jl new file mode 100644 index 0000000..a35506e --- /dev/null +++ b/src/algorithms/diamondsquare.jl @@ -0,0 +1,206 @@ +""" + DiamondSquare <: NeutralLandscapeMaker + + DiamondSquare(; H = 0.5) + DiamondSquare(H) + +This type generates a neutral landscape using the diamond-squares +algorithm, which produces fractals with variable spatial autocorrelation. + +https://en.wikipedia.org/wiki/Diamond-square_algorithm + +The algorithm is named diamond-square because it is an iterative procedure of +"diamond" and "square" steps. + + +The degree of spatial autocorrelation is controlled by a parameter `H`, +which varies from 0.0 (low autocorrelation) to 1.0 (high autocorrelation) --- +note this is non-inclusive and H = 0 and H = 1 will not behave as expected. +The result of the diamond-square algorithm is a fractal with dimension D = 2 + H. + +A similar algorithm, midpoint-displacement, functions almost +identically, except that in DiamondSquare, the square step interpolates +edge midpoints from the nearest two corners and the square's center, where as +midpoint-displacement only interpolates from the nearest corners (see `MidpointDisplacement`). + +""" +@kwdef struct DiamondSquare <: NeutralLandscapeMaker + H::Float64 = 0.5 + function DiamondSquare(H::T) where {T <: Real} + @assert 0 <= H < 1 + new(H) + end +end + +""" + MidpointDisplacement <: NeutralLandscapeMaker + + MidpointDisplacement(; H = 0.5) + +Creates a midpoint-displacement algorithm object `MidpointDisplacement`. +The degree of spatial autocorrelation is controlled by a parameter `H`, +which varies from 0.0 (low autocorrelation) to 1.0 (high autocorrelation) --- +note this is non-inclusive and H = 0 and H = 1 will not behave as expected. + +A similar algorithm, diamond-square, functions almost +identically, except that in diamond-square, the square step interpolates +edge midpoints from the nearest two corners and the square's center, where as +`MidpointDisplacement` only interpolates from the nearest corners (see `DiamondSquare`). +""" +@kwdef struct MidpointDisplacement <: NeutralLandscapeMaker + H::Float64 = 0.5 + function MidpointDisplacement(H::T) where {T <: Real} + @assert 0 <= H < 1 + new(H) + end +end + +# Check if `mat` is the right size. If mat is not the correct size (DiamondSquare +# can only run on a lattice of size NxN where N = (2^n)+1 for integer n), +# allocates the smallest lattice large enough to contain `mat` that can run +# DiamondSquare. +function _landscape!(mat, alg::Union{DiamondSquare, MidpointDisplacement}; kw...) where {IT <: Integer} + + rightSize::Bool = _isPowerOfTwo(size(mat)[1]-1) && _isPowerOfTwo(size(mat)[2]-1) + latticeSize::Int = size(mat)[1] + + dsMat = mat + if !rightSize + dim1, dim2 = size(mat) + smallestContainingLattice::Int = 2^ceil(log2(max(dim1, dim2))) + 1 + dsMat = zeros(smallestContainingLattice, smallestContainingLattice) + end + _diamondsquare!(dsMat, alg) + + mat .= dsMat[1:size(mat)[1], 1:size(mat)[2]] +end + +# Runs the diamond-square algorithm on a matrix `mat` of size +# `NxN`, where `N=(2^n)+1` for some integer `n`, i.e (N=5,9,17,33,65) + +# Diamond-square is an iterative procedure, where the lattice is divided +# into subsquares in subsequent rounds. At each round, the subsquares shrink in size, +# as previously uninitialized values in the lattice are interpolated as a mean of nearby points plus random displacement. +# As the rounds increase, the magnitude of this displacement decreases. This creates spatioautocorrelation, which is controlled +# by a single parameter `H` which varies between `0` (no autocorrelation) and `1` (high autocorrelation) +function _diamondsquare!(mat, alg) + latticeSize = size(mat)[1] + numberOfRounds::Int = log2(latticeSize-1) + _initializeDiamondSquare!(mat, alg) + + for round in 0:(numberOfRounds-1) # counting from 0 saves us a headache later + subsquareSideLength::Int = 2^(numberOfRounds-(round)) + numberOfSubsquaresPerAxis::Int = ((latticeSize-1) / subsquareSideLength)-1 + + for x in 0:numberOfSubsquaresPerAxis # iterate over the subsquares within the lattice at this side length + for y in 0:numberOfSubsquaresPerAxis + subsquareCorners = _subsquareCornerCoordinates(x,y,subsquareSideLength) + + _diamond!(mat, alg, round, subsquareCorners) + _square!(mat, alg, round, subsquareCorners) + end + end + end +end + +# Initialize's the `DiamondSquare` algorithm by displacing the four corners of the +# lattice using `displace`, scaled by the algorithm's autocorrelation `H`. +function _initializeDiamondSquare!(mat, alg) + latticeSize = size(mat)[1] + corners = _subsquareCornerCoordinates(0,0, latticeSize-1) + for mp in corners + mat[mp...] = _displace(alg.H, 1) + end +end + +# Returns the coordinates for the corners of the subsquare (x,y) given a side-length `sideLength`. +function _subsquareCornerCoordinates(x::Int, y::Int, sideLength::Int) + corners = [1 .+ sideLength.*i for i in [(x,y), (x+1, y), (x, y+1), (x+1, y+1)]] +end + +# Runs the diamond step of the `DiamondSquare` algorithm on the square defined by +# `corners` on the matrix `mat`. The center of the square is interpolated from the +# four corners, and is displaced. The displacement is drawn according to `alg.H` and round using `displace` +function _diamond!(mat, alg, round::Int, corners::AbstractVector{Tuple{Int, Int}}) + centerPt = _centerCoordinate(corners) + mat[centerPt...] = _interpolate(mat, corners) + _displace(alg.H, round) +end + +# Runs the square step of the `DiamondSquare` algorithm on the square defined +# by `corners` on the matrix `mat`. The midpoint of each edge of this square is interpolated +# by computing the mean value of the two corners on the edge and the center of the square, and the +# displacing it. The displacement is drawn according to `alg.H` and round using `displace` +function _square!(mat, alg::DiamondSquare, round::Int, corners::AbstractVector{Tuple{Int, Int}}) + bottomLeft,bottomRight,topLeft,topRight = corners + leftEdge, bottomEdge, topEdge, rightEdge = _edgeMidpointCoordinates(corners) + centerPoint = _centerCoordinate(corners) + + mat[leftEdge...] = _interpolate(mat, [topLeft,bottomLeft,centerPoint]) + _displace(alg.H, round) + mat[bottomEdge...] = _interpolate(mat, [bottomLeft,bottomRight,centerPoint]) + _displace(alg.H, round) + mat[topEdge...] = _interpolate(mat, [topLeft,topRight,centerPoint]) + _displace(alg.H, round) + mat[rightEdge...] = _interpolate(mat, [topRight,bottomRight,centerPoint]) + _displace(alg.H, round) +end + +# Runs the square step of the `MidpointDisplacement` algorithm on the square defined +# by `corners` on the matrix `mat`. The midpoint of each edge of this square is interpolated +# by computing the mean value of the two corners on the edge and the center of the square, and the +# displacing it. The displacement is drawn according to `alg.H` and round using `displace` +function _square!(mat, alg::MidpointDisplacement, round::Int, corners::AbstractVector{Tuple{Int, Int}}) + bottomLeft,bottomRight,topLeft,topRight = corners + leftEdge, bottomEdge, topEdge, rightEdge = _edgeMidpointCoordinates(corners) + mat[leftEdge...] = _interpolate(mat, [topLeft,bottomLeft]) + _displace(alg.H, round) + mat[bottomEdge...] = _interpolate(mat, [bottomLeft,bottomRight]) + _displace(alg.H, round) + mat[topEdge...] = _interpolate(mat, [topLeft,topRight]) + _displace(alg.H, round) + mat[rightEdge...] = _interpolate(mat, [topRight,bottomRight]) + _displace(alg.H, round) +end + +# Computes the mean of a set of points, represented as a list of indicies to a matrix `mat`. +function _interpolate(mat, points::AbstractVector{Tuple{Int,Int}}) + return mean(mat[pt...] for pt in points) +end + +# `displace` produces a random value as a function of `H`, which is the +# autocorrelation parameter used in `DiamondSquare` and must be between `0` +# and `1`, and `round` which describes the current tiling size for the +# DiamondSquare() algorithm. + +# Random value are drawn from a Gaussian distribution using `Distribution.Normal` +# The standard deviation of this Gaussian, σ, is set to (1/2)^(round*H), which will +# move from 1.0 to 0 as `round` increases. +function _displace(H::Float64, round::Int) + σ = 0.5^(round*H) + return(rand(Normal(0, σ))) +end + +# Returns the center coordinate for a square defined by `corners` for the +# `DiamondSquare` algorithm. +function _centerCoordinate(corners::AbstractVector{Tuple{Int,Int}}) + bottomLeft,bottomRight,topLeft,topRight = corners + centerX::Int = (_xcoord(bottomLeft)+_xcoord(bottomRight)) ÷ 2 + centerY::Int = (_ycoord(topRight)+_ycoord(bottomRight)) ÷ 2 + return (centerX, centerY) +end + +# Returns the x-coordinate from a lattice coordinate `pt`. +_xcoord(pt::Tuple{Int,Int}) = pt[1] +# Returns the y-coordinate from a lattice coordinate `pt`. +_ycoord(pt::Tuple{Int,Int}) = pt[2] + +# Returns an array of midpoints for a square defined by `corners` for the `DiamondSquare` algorithm. +function _edgeMidpointCoordinates(corners::AbstractVector{Tuple{Int,Int}}) + # bottom left, bottom right, top left, top right + bottomLeft,bottomRight,topLeft,topRight = corners + + leftEdgeMidpoint::Tuple{Int,Int} = (_xcoord(bottomLeft), (_ycoord(bottomLeft)+_ycoord(topLeft))÷2 ) + bottomEdgeMidpoint::Tuple{Int,Int} = ( (_xcoord(bottomLeft)+ _xcoord(bottomRight))÷2, _ycoord(bottomLeft) ) + topEdgeMidpoint::Tuple{Int,Int} = ( (_xcoord(topLeft)+_xcoord(topRight))÷2, _ycoord(topLeft)) + rightEdgeMidpoint::Tuple{Int,Int} = ( _xcoord(bottomRight), (_ycoord(bottomRight)+_ycoord(topRight))÷2) + + edgeMidpoints = [leftEdgeMidpoint, bottomEdgeMidpoint, topEdgeMidpoint, rightEdgeMidpoint] + return edgeMidpoints +end + +# Determines if `x`, an integer, can be expressed as `2^n`, where `n` is also an integer. +function _isPowerOfTwo(x::IT) where {IT <: Integer} + return x & (x-1) == 0 +end diff --git a/src/algorithms/discretevoronoi.jl b/src/algorithms/discretevoronoi.jl new file mode 100644 index 0000000..8155083 --- /dev/null +++ b/src/algorithms/discretevoronoi.jl @@ -0,0 +1,21 @@ +""" + DiscreteVoronoi <: NeutralLandscapeMaker + + DiscreteVoronoi(; n=3) + DiscreteVoronoi(n) + +This type provides a rasterization of a Voronoi-like diagram. +Assigns a value to each patch using a 1-NN algorithmm with `n` initial clusters. +It is a `NearestNeighborElement` algorithmm with `k` neighbors set to 1. +The default is to use three clusters. +""" +@kwdef struct DiscreteVoronoi <: NeutralLandscapeMaker + n::Int64 = 3 + function DiscreteVoronoi(n::Int64) + @assert n > 0 + new(n) + end +end + +_landscape!(mat, alg::DiscreteVoronoi) = _landscape!(mat, NearestNeighborElement(alg.n, 1)) + diff --git a/src/algorithms/distancegradient.jl b/src/algorithms/distancegradient.jl new file mode 100644 index 0000000..d2e2b1d --- /dev/null +++ b/src/algorithms/distancegradient.jl @@ -0,0 +1,24 @@ +""" + DistanceGradient <: NeutralLandscapeMaker + + DistanceGradient(; sources=[1]) + DistanceGradient(sources) + +The `sources` field is a `Vector{Integer}` of *linear* indices of the matrix, +from which the distance must be calculated. +""" +@kwdef struct DistanceGradient <: NeutralLandscapeMaker + sources::Vector{Integer} = [1] +end + +function _landscape!(mat, alg::DistanceGradient) + @assert maximum(alg.sources) <= length(mat) + @assert minimum(alg.sources) > 0 + coordinates = Matrix{Float64}(undef, (2, prod(size(mat)))) + for (i, p) in enumerate(Iterators.product(axes(mat)...)) + coordinates[1:2, i] .= p + end + tree = KDTree(coordinates[:,alg.sources]) + mat[:] .= nn(tree, coordinates)[2] + return mat +end diff --git a/src/algorithms/edgegradient.jl b/src/algorithms/edgegradient.jl new file mode 100644 index 0000000..e0839f5 --- /dev/null +++ b/src/algorithms/edgegradient.jl @@ -0,0 +1,21 @@ +""" + EdgeGradient <: NeutralLandscapeMaker + + EdgeGradient(; direction=360rand()) + EdgeGradient(direction) + +This type is used to generate an edge gradient landscape, where values change +as a bilinear function of the *x* and *y* coordinates. The direction is +expressed as a floating point value, which will be in *[0,360]*. The inner +constructor takes the mod of the value passed and 360, so that a value that is +out of the correct interval will be corrected. +""" +@kwdef struct EdgeGradient <: NeutralLandscapeMaker + direction::Float64 = 360rand() + EdgeGradient(x::T) where {T <: Real} = new(mod(x, 360.0)) +end + +function _landscape!(mat, alg::EdgeGradient) where {IT <: Integer} + _landscape!(mat, PlanarGradient(alg.direction)) + mat .= -2.0abs.(0.5 .- mat) .+ 1.0 +end diff --git a/src/algorithms/nncluster.jl b/src/algorithms/nncluster.jl new file mode 100644 index 0000000..a52af29 --- /dev/null +++ b/src/algorithms/nncluster.jl @@ -0,0 +1,33 @@ +""" + NearestNeighborCluster <: NeutralLandscapeMaker + + NearestNeighborCluster(; p=0.5, n=:rook) + NearestNeighborCluster(p, [n=:rook]) + +Create a random cluster nearest-neighbour neutral landscape model with +values ranging 0-1. `p` sets the density of original clusters, and `n` +sets the neighborhood for clustering (see `?label` for neighborhood options) +""" +@kwdef struct NearestNeighborCluster <: NeutralLandscapeMaker + p::Float64 = 0.5 + n::Symbol = :rook + function NearestNeighborCluster(p::Float64, n::Symbol = :rook) + @assert p > 0 + @assert n ∈ (:rook, :queen, :diagonal) + new(p,n) + end +end + +function _landscape!(mat, alg::NearestNeighborCluster) + _landscape!(mat, NoGradient()) + classify!(mat, [alg.p, 1 - alg.p]) + replace!(mat, 2.0 => NaN) + clusters, nClusters = label(mat, alg.n) + coordinates = _coordinatematrix(clusters) + sources = findall(!isnan, vec(clusters)) + tree = KDTree(coordinates[:,sources]) + clusters[:] .= clusters[sources[nn(tree, coordinates)[1]]] + randvals = rand(nClusters) + mat .= randvals[Int.(clusters)] + mat +end diff --git a/src/algorithms/nnelement.jl b/src/algorithms/nnelement.jl new file mode 100644 index 0000000..a5ba8cf --- /dev/null +++ b/src/algorithms/nnelement.jl @@ -0,0 +1,43 @@ +""" + NearestNeighborElement <: NeutralLandscapeMaker + + NearestNeighborElement(; n=3, k=1) + NearestNeighborElement(n, [k=1]) + +Assigns a value to each patch using a k-NN algorithmm with `n` initial clusters +and `k` neighbors. The default is to use three cluster and a single neighbor. +""" +@kwdef struct NearestNeighborElement <: NeutralLandscapeMaker + n::Int64 = 3 + k::Int64 = 1 + function NearestNeighborElement(n::Int64, k::Int64 = 1) + @assert n > 0 + @assert k > 0 + @assert k <= n + new(n,k) + end +end + +function _coordinatematrix(mat) + coordinates = Matrix{Float64}(undef, (2, prod(size(mat)))) + for (i, p) in enumerate(Iterators.product(axes(mat)...)) + coordinates[1:2, i] .= p + end + coordinates +end + +function _landscape!(mat, alg::NearestNeighborElement) + fill!(mat, 0.0) + clusters = sample(eachindex(mat), alg.n; replace=false) + for (i,n) in enumerate(clusters) + mat[n] = i + end + coordinates = _coordinatematrix(mat) + tree = KDTree(coordinates[:,clusters]) + if alg.k == 1 + mat[:] .= nn(tree, coordinates)[1] + else + mat[:] .= map(mean, knn(tree, coordinates, alg.k)[1]) + end + return mat +end diff --git a/src/algorithms/nogradient.jl b/src/algorithms/nogradient.jl new file mode 100644 index 0000000..c05f396 --- /dev/null +++ b/src/algorithms/nogradient.jl @@ -0,0 +1,11 @@ +""" + NoGradient <: NeutralLandscapeMaker + + NoGradient() + +This type is used to generate a random landscape with no gradients +""" +struct NoGradient <: NeutralLandscapeMaker +end + +_landscape!(mat, ::NoGradient) = rand!(mat) diff --git a/src/algorithms/perlinnoise.jl b/src/algorithms/perlinnoise.jl new file mode 100644 index 0000000..3265682 --- /dev/null +++ b/src/algorithms/perlinnoise.jl @@ -0,0 +1,136 @@ +""" + PerlinNoise <: NeutralLandscapeMaker + + PerlinNoise(; kw...) + PerlinNoise(periods, [octaves=1, lacunarity=2, persistance=0.5, valley=:u]) + +Create a Perlin noise neutral landscape model with values ranging 0-1. + +# Keywords + +- `periods::Tuple{Int,Int}=(1,1)`: the number of periods of Perlin noise across row and + column dimensions for the first octave. +- `octaves::Int=1`: the number of octaves that will form the Perlin noise. +- `lacunarity::Int=2` : the rate at which the frequency of periods increases for each + octive. +- `persistance::Float64=0.5` : the rate at which the amplitude of periods decreases for each + octive. +- `valley::Symbol=`:u`: the kind of valley bottom that will be mimicked: `:u` produces + u-shaped valleys, `:v` produces v-shaped valleys, and `:-` produces flat bottomed + valleys. + +Note: This is a memory-intensive algorithm with some settings. Be careful using larger +prime numbers for `period` when also using a large array size, high lacuarity and/or many +octaves. Memory use scales with the lowest common multiple of `periods`. +""" +@kwdef struct PerlinNoise <: NeutralLandscapeMaker + periods::Tuple{Int,Int} = (1, 1) + octaves::Int = 1 + lacunarity::Int = 2 + persistance::Float64 = 0.5 + valley::Symbol = :u + function PerlinNoise(periods, octaves=1, lacunarity=2, persistance=0.5, valley=:u) + new(periods, octaves, lacunarity, persistance, valley) + end +end + +function _landscape!(A, alg::PerlinNoise) + # nRow must equal nCol so determine the dimension of the smallest square + dim = max(size(A)...) + # Check the dim is a multiple of each octives maximum number of periods and + # expand dim if needed + rperiodsmax = alg.periods[1] * alg.lacunarity^(alg.octaves - 1) + cperiodsmax = alg.periods[2] * alg.lacunarity^(alg.octaves - 1) + periodsmultiple = lcm(rperiodsmax, cperiodsmax) # lowest common multiple + if dim % periodsmultiple != 0 + dim = ceil(Int, dim / periodsmultiple) * periodsmultiple + end + + # Generate the Perlin noise + noise = zeros(eltype(A), (dim, dim)) + meshbuf = Array{eltype(A),3}(undef, dim, dim, 2) + nbufs = ntuple(_->Array{eltype(A),2}(undef, dim, dim), 4) + for octave in 0:alg.octaves-1 + octave_noise!(noise, meshbuf, nbufs, alg, octave, dim, dim) + end + + # Randomly extract the desired array size + noiseview = _view_from_square(noise, size(A)...) + + # Rescale the Perlin noise to mimic different kinds of valley bottoms + return if alg.valley == :u + A .= noiseview + elseif alg.valley == :v + A .= abs.(noiseview) + elseif alg.valley == :- + A .= noiseview.^2 + else + error("$(alg.valley) not recognised for `valley` use `:u`, `:v` or `:-`") + end +end + +function octave_noise!( + noise, mesh, (n11, n21, n12, n22), alg::PerlinNoise, octave, nrow, ncol +) + f(t) = @fastmath 6 * t ^ 5 - 15 * t ^ 4 + 10 * t ^ 3 # Wut + + # Mesh + rp, cp = alg.periods .* alg.lacunarity^(octave) + delta = (rp / nrow, cp / ncol) + ranges = range(0, rp-delta[1], length=nrow), range(0, cp-delta[2], length=ncol) + _mesh!(mesh, ranges...) + @fastmath mesh .%= 1 + + # Gradients + # This allocates, but the gradients size changes with octave so it needs to + # be a very smart in-place `repeat!` or some kind of generator, and the improvement + # may not be than large (~20%). + angles = 2pi .* rand(rp + 1, cp + 1) + @fastmath gradients = cat(cos.(angles), sin.(angles); dims=3) + d = (nrow ÷ rp, ncol ÷ cp) + grad = repeat(gradients, inner=[d[1], d[2], 1]) + + g111 = @view grad[1:ncol, 1:ncol, 1] + g211 = @view grad[end-nrow+1:end, 1:ncol, 1] + g121 = @view grad[1:ncol, end-ncol+1:end, 1] + g221 = @view grad[end-nrow+1:end, end-ncol+1:end, 1] + g112 = @view grad[1:ncol, 1:ncol, 2] + g212 = @view grad[end-nrow+1:end, 1:ncol, 2] + g122 = @view grad[1:ncol, end-ncol+1:end, 2] + g222 = @view grad[end-nrow+1:end, end-ncol+1:end, 2] + + # Ramps + m1 = @view mesh[:, :, 1] + m2 = @view mesh[:, :, 2] + n11 .= ((m1 .+ m2 ) .* g111 .+ (m1 .+ m2 ) .* g112) + n21 .= ((m1 .-1 .+ m2 ) .* g211 .+ (m1 .-1 .+ m2 ) .* g212) + n12 .= ((m1 .+ m2 .- 1) .* g121 .+ (m1 .+ m2 .- 1) .* g122) + n22 .= ((m1 .-1 .+ m2 .- 1) .* g221 .+ (m1 .-1 .+ m2 .- 1) .* g222) + + # Interpolation + mesh .= f.(mesh) + t1 = @view mesh[:, :, 1] + t2 = @view mesh[:, :, 2] + noise .+= sqrt(2) .* (alg.persistance ^ octave) .* + ((1 .- t2) .* (n11 .* (1 .- t1) .+ t1 .* n21) .+ + t2 .* (n12 .* (1 .- t1) .+ t1 .* n22)) + return noise +end + +function _mesh!(A, x, y) + for (i, ival) in enumerate(x), j in 1:length(y) + A[i, j, 1] = ival + end + for i in 1:length(x), (j, jval) in enumerate(y) + A[i, j, 2] = jval + end + return A +end + +function _view_from_square(source, nrow, ncol) + # Extract a portion of the array to match the dimensions + dim = size(source, 1) + startrow = rand(1:(dim - nrow + 1)) + startcol = rand(1:(dim - ncol + 1)) + return @view source[startrow:startrow + nrow - 1, startcol:startcol + ncol - 1] +end diff --git a/src/algorithms/planargradient.jl b/src/algorithms/planargradient.jl new file mode 100644 index 0000000..601b06b --- /dev/null +++ b/src/algorithms/planargradient.jl @@ -0,0 +1,24 @@ +""" + PlanarGradient <: NeutralLandscapeMaker + + PlanarGradient(; direction=360rand()) + PlanarGradient(direction) + +This type is used to generate a planar gradient landscape, where values change +as a bilinear function of the *x* and *y* coordinates. The direction is +expressed as a floating point value, which will be in *[0,360]*. The inner +constructor takes the mod of the value passed and 360, so that a value that is +out of the correct interval will be corrected. +""" +@kwdef struct PlanarGradient <: NeutralLandscapeMaker + direction::Float64 = 360rand() + PlanarGradient(x::T) where {T <: Real} = new(mod(x, 360.0)) +end + +function _landscape!(mat, alg::PlanarGradient) where {IT <: Integer} + eastness = sin(deg2rad(alg.direction)) + southness = -1cos(deg2rad(alg.direction)) + rows, cols = axes(mat) + mat .= southness .* rows .+ (eastness .* cols)' + return mat +end diff --git a/src/algorithms/rectangularcluster.jl b/src/algorithms/rectangularcluster.jl new file mode 100644 index 0000000..f2b6ff6 --- /dev/null +++ b/src/algorithms/rectangularcluster.jl @@ -0,0 +1,29 @@ +""" + RectangularCluster <: NeutralLandscapeMaker + + RectangularCluster(; minimum=2, maximum=4) + RectangularCluster(minimum, [maximum=4]) + +Fills the landscape with rectangles containing a random value. The size of each +rectangle/patch is between `minimum` and `maximum` (the two can be equal for a +fixed size rectangle). +""" +@kwdef struct RectangularCluster <: NeutralLandscapeMaker + minimum::Integer = 2 + maximum::Integer = 4 + function RectangularCluster(minimum::T, maximum::T=4) where {T <: Integer} + @assert 0 < minimum <= maximum + new(minimum, maximum) + end +end + +function _landscape!(mat, alg::RectangularCluster) + mat .= -1.0 + while minimum(mat) == -1.0 + width, height = rand(alg.minimum:alg.maximum, 2) + row = rand(1:(size(mat,1)-(width-1))) + col = rand(1:(size(mat,2)-(height-1))) + mat[row:(row+(width-1)) , col:(col+(height-1))] .= rand() + end + return mat +end diff --git a/src/algorithms/wavesurface.jl b/src/algorithms/wavesurface.jl new file mode 100644 index 0000000..fa86d82 --- /dev/null +++ b/src/algorithms/wavesurface.jl @@ -0,0 +1,22 @@ +""" + WaveSurface <: NeutralLandscapeMaker + + WaveSurface(; direction=360rand(), periods=1) + WaveSurface(direction, [periods=1]) + +Creates a sinusoidal landscape with a `direction` and a number of `periods`. If +neither are specified, there will be a single period of random direction. +""" +@kwdef struct WaveSurface <: NeutralLandscapeMaker + direction::Float64 = 360rand() + periods::Int64 = 1 + function WaveSurface(direction::T, periods::K = 1) where {T <: Real, K <: Integer} + @assert periods >= 1 + new(mod(direction, 360.0), periods) + end +end + +function _landscape!(mat, alg::WaveSurface) + rand!(mat, PlanarGradient(alg.direction)) + mat .= sin.(mat .* (2π * alg.periods)) +end From c8afa6d8c91707f66d0bddb823e443354fe16b91 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Fri, 31 Dec 2021 16:51:28 -0500 Subject: [PATCH 07/31] splint rate into rate and direction --- src/updaters/update.jl | 28 +++++++++++++++------------- 1 file changed, 15 insertions(+), 13 deletions(-) diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 3b4a51d..f582726 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -15,33 +15,35 @@ abstract type NeutralLandscapeUpdater end which defines the expected change in any cell per timestep. """ rate(up::NeutralLandscapeUpdater) = up.rate - -generator(up::NeutralLandscapeUpdater) = up.generator +spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater +direction(up::NeutralLandscapeUpdater) = up.direction function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) end -@kwdef struct SpatiallyAutocorrelatedUpdater{G,R} <: NeutralLandscapeUpdater - generator::G = DiamondSquare(0.5) - rate::R = 0.1 +@kwdef struct SpatiallyAutocorrelatedUpdater{SU,D,S} <: NeutralLandscapeUpdater + spatialupdater::SU = DiamondSquare(0.5) + direction::D = 0.1 + rate::S = 0.1 end function _update(sau::SpatiallyAutocorrelatedUpdater, mat) - change = rand(generator(sau), size(mat)) - delta = rate(sau) .+ rate(sau) .* transform(fit(ZScoreTransform, change), change) + change = rand(spatialupdater(sau), size(mat)) + delta = direction(sau) .+ rate(sau) .* transform(fit(ZScoreTransform, change), change) mat .+ delta end -@kwdef struct SpatiotemporallyAutocorrelatedUpdater{G,R,T} <: NeutralLandscapeUpdater - generator::G = DiamondSquare(0.1) - rate::R = 0.1 - temporalupdater::T = BrownianMotion() +@kwdef struct SpatiotemporallyAutocorrelatedUpdater{SU,TU,D,R} <: NeutralLandscapeUpdater + spatialupdater::SU = DiamondSquare(0.1) + temporalupdater::TU = BrownianMotion() + direction::D = 0.1 + rate::R = 0.1 end function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) - change = rand(generator(stau), size(mat)) + change = rand(spatialupdater(stau), size(mat)) temporalshift = broadcast(x->update(stau.temporalupdater, x), mat) z = transform(fit(ZScoreTransform, change), change) - delta = rate(stau) .+ (rate(stau) .* z) .+ (temporalshift .* z) + delta = direction(stau) .+ (rate(stau) .* z) .+ (temporalshift .* z) mat .+ delta end From 86c67b18d7ce8c9cdbba0dee204cfe240f2b142f Mon Sep 17 00:00:00 2001 From: michael catchen Date: Sat, 8 Jan 2022 12:05:05 -0500 Subject: [PATCH 08/31] updates --- src/NeutralLandscapes.jl | 2 +- src/updaters/update.jl | 4 ++++ 2 files changed, 5 insertions(+), 1 deletion(-) diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index a3c6274..e896a90 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -43,7 +43,7 @@ include("algorithms/wavesurface.jl") include("updaters/update.jl") include("updaters/temporal.jl") export TemporalUpdater, BrownianMotion, OhrsteinUhlenbeck -export NeutralLandscapeUpdater, update, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater +export NeutralLandscapeUpdater, update, update!, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater end # module diff --git a/src/updaters/update.jl b/src/updaters/update.jl index f582726..ef8d9ce 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -21,6 +21,10 @@ direction(up::NeutralLandscapeUpdater) = up.direction function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) end +function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} + mat = _update(updater, mat) +end + @kwdef struct SpatiallyAutocorrelatedUpdater{SU,D,S} <: NeutralLandscapeUpdater spatialupdater::SU = DiamondSquare(0.5) direction::D = 0.1 From 9dd2f7372ecb54f55352b76a06e3f3fdde2bedb8 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Sat, 8 Jan 2022 13:49:37 -0500 Subject: [PATCH 09/31] :wrench: nearly there --- src/NeutralLandscapes.jl | 5 +++-- src/updaters/spatial.jl | 12 +++++++++++ src/updaters/spatiotemporal.jl | 14 +++++++++++++ src/updaters/temporal.jl | 24 +++++++-------------- src/updaters/update.jl | 38 ---------------------------------- 5 files changed, 37 insertions(+), 56 deletions(-) create mode 100644 src/updaters/spatial.jl create mode 100644 src/updaters/spatiotemporal.jl diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index e896a90..89139d6 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -42,8 +42,9 @@ include("algorithms/wavesurface.jl") include("updaters/update.jl") include("updaters/temporal.jl") -export TemporalUpdater, BrownianMotion, OhrsteinUhlenbeck -export NeutralLandscapeUpdater, update, update!, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater +include("updaters/spatial.jl") +include("updaters/spatiotemporal.jl") +export update, update!, NeutralLandscapeUpdater, TemporallyAutocorrelatedUpdater, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater end # module diff --git a/src/updaters/spatial.jl b/src/updaters/spatial.jl new file mode 100644 index 0000000..695f3e0 --- /dev/null +++ b/src/updaters/spatial.jl @@ -0,0 +1,12 @@ + + +@kwdef struct SpatiallyAutocorrelatedUpdater{SU,D,S} <: NeutralLandscapeUpdater + spatialupdater::SU = DiamondSquare(0.5) + direction::D = 0.1 + rate::S = 0.1 +end +function _update(sau::SpatiallyAutocorrelatedUpdater, mat) + change = rand(spatialupdater(sau), size(mat)) + delta = direction(sau) .+ rate(sau) .* transform(fit(ZScoreTransform, change), change) + mat .+ delta +end diff --git a/src/updaters/spatiotemporal.jl b/src/updaters/spatiotemporal.jl new file mode 100644 index 0000000..b657301 --- /dev/null +++ b/src/updaters/spatiotemporal.jl @@ -0,0 +1,14 @@ + + +@kwdef struct SpatiotemporallyAutocorrelatedUpdater{SU,D,R} <: NeutralLandscapeUpdater + spatialupdater::SU = DiamondSquare(0.1) + direction::D = 0.1 + rate::R = 0.1 +end +function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) + change = rand(spatialupdater(stau), size(mat)) + temporalshift = rand(Normal(), size(mat)) + z = transform(fit(ZScoreTransform, change), change) + delta = direction(stau) .+ (rate(stau) .* z) .+ (rate(stau) .* temporalshift) + mat .+ delta +end \ No newline at end of file diff --git a/src/updaters/temporal.jl b/src/updaters/temporal.jl index cbf296c..89a3861 100644 --- a/src/updaters/temporal.jl +++ b/src/updaters/temporal.jl @@ -1,19 +1,11 @@ -abstract type TemporalUpdater end -update(tu::T, x) where {T<:TemporalUpdater} = _getchange(tu, x) - -@kwdef struct BrownianMotion{S} <: TemporalUpdater - σ::S = 0.1 +@kwdef struct TemporallyAutocorrelatedUpdater{D,S} <: NeutralLandscapeUpdater + spatialupdater = missing + direction::D = 0.1 + rate::S = 0.1 end - -_getchange(bm::BrownianMotion, x) = rand(Normal(bm.σ)) - - -@kwdef struct OhrsteinUhlenbeck{S,T,F} <: TemporalUpdater - θ::T = 1. - σ::S = 0.1 - dt::F = 0.01 +function _update(sau::TemporallyAutocorrelatedUpdater, mat) + U = rand(Normal(), size(mat)) + Δ = direction(sau) .+ rate(sau) .* U + mat .+ Δ end - -_getchange(ou::OhrsteinUhlenbeck, x) = ou.dt*(-ou.θ*x + rand(Normal(ou.σ))) - diff --git a/src/updaters/update.jl b/src/updaters/update.jl index ef8d9ce..3ab338a 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -25,41 +25,3 @@ function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} mat = _update(updater, mat) end -@kwdef struct SpatiallyAutocorrelatedUpdater{SU,D,S} <: NeutralLandscapeUpdater - spatialupdater::SU = DiamondSquare(0.5) - direction::D = 0.1 - rate::S = 0.1 -end -function _update(sau::SpatiallyAutocorrelatedUpdater, mat) - change = rand(spatialupdater(sau), size(mat)) - delta = direction(sau) .+ rate(sau) .* transform(fit(ZScoreTransform, change), change) - mat .+ delta -end - - -@kwdef struct SpatiotemporallyAutocorrelatedUpdater{SU,TU,D,R} <: NeutralLandscapeUpdater - spatialupdater::SU = DiamondSquare(0.1) - temporalupdater::TU = BrownianMotion() - direction::D = 0.1 - rate::R = 0.1 -end -function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) - change = rand(spatialupdater(stau), size(mat)) - temporalshift = broadcast(x->update(stau.temporalupdater, x), mat) - z = transform(fit(ZScoreTransform, change), change) - delta = direction(stau) .+ (rate(stau) .* z) .+ (temporalshift .* z) - mat .+ delta -end - - - - -""" - example of how its used - -""" -#env = rand(MidpointDisplacement(0.9), (250,250)) -#sau = SpatiallyAutocorrelatedUpdater() -#new = update!(sau, env) -#plot(heatmap.([env, new])...) - From a39fbec0a9a936d139a991b244d1f613cb2f272f Mon Sep 17 00:00:00 2001 From: michael catchen Date: Mon, 10 Jan 2022 11:23:37 -0500 Subject: [PATCH 10/31] up --- Project.toml | 1 + src/updaters/update.jl | 11 ++++++++++- 2 files changed, 11 insertions(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index e8224ae..f79b0ed 100644 --- a/Project.toml +++ b/Project.toml @@ -8,6 +8,7 @@ DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" NearestNeighbors = "b8a86587-4115-5ab1-83bc-aa920d37bbce" +ProfileView = "c46f51b8-102a-5cf2-8d2c-8597cb0e0da7" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 3ab338a..0e9f4ca 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -18,10 +18,19 @@ rate(up::NeutralLandscapeUpdater) = up.rate spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater direction(up::NeutralLandscapeUpdater) = up.direction + function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) end +function update(updater::T, mat, n::I) where {T<:NeutralLandscapeUpdater, I<:Integer} + sequence = [zeros(size(mat)) for _ in 1:n] + sequence[begin] = mat + for i in 2:n + sequence[i] = _update(updater, sequence[i-1]) + end + sequence +end function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} - mat = _update(updater, mat) + mat[:] = _update(updater, mat) end From 1581042f890b749bde24e883eb7650a8ead61af6 Mon Sep 17 00:00:00 2001 From: rafaqz Date: Mon, 10 Jan 2022 12:20:05 +0100 Subject: [PATCH 11/31] performance improvements --- src/NeutralLandscapes.jl | 2 +- src/algorithms/diamondsquare.jl | 34 ++++++++++++------------- src/algorithms/distancegradient.jl | 23 +++++++++++++---- src/algorithms/nncluster.jl | 20 +++++++++++---- src/algorithms/nnelement.jl | 37 +++++++++++++++------------- src/algorithms/perlinnoise.jl | 31 +++++++++++------------ src/algorithms/planargradient.jl | 2 +- src/algorithms/rectangularcluster.jl | 2 +- src/classify.jl | 6 ++--- src/landscape.jl | 7 +++--- test/rand.jl | 4 +++ 11 files changed, 97 insertions(+), 71 deletions(-) diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index bb8a7d6..81d0904 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -5,7 +5,7 @@ using StatsBase: sample using Random: rand! using Statistics: quantile, mean using Distributions: Normal -using NearestNeighbors: KDTree, knn, nn +using NearestNeighbors: KDTree, knn, nn, always_false, knn_point!, SVector using DataStructures: IntDisjointSets, union!, find_root, push! using Base: @kwdef diff --git a/src/algorithms/diamondsquare.jl b/src/algorithms/diamondsquare.jl index a35506e..a9c7345 100644 --- a/src/algorithms/diamondsquare.jl +++ b/src/algorithms/diamondsquare.jl @@ -115,13 +115,13 @@ end # Returns the coordinates for the corners of the subsquare (x,y) given a side-length `sideLength`. function _subsquareCornerCoordinates(x::Int, y::Int, sideLength::Int) - corners = [1 .+ sideLength.*i for i in [(x,y), (x+1, y), (x, y+1), (x+1, y+1)]] + return (1 .+ sideLength.*i for i in ((x,y), (x+1, y), (x, y+1), (x+1, y+1))) end # Runs the diamond step of the `DiamondSquare` algorithm on the square defined by # `corners` on the matrix `mat`. The center of the square is interpolated from the # four corners, and is displaced. The displacement is drawn according to `alg.H` and round using `displace` -function _diamond!(mat, alg, round::Int, corners::AbstractVector{Tuple{Int, Int}}) +function _diamond!(mat, alg, round::Int, corners) centerPt = _centerCoordinate(corners) mat[centerPt...] = _interpolate(mat, corners) + _displace(alg.H, round) end @@ -130,32 +130,32 @@ end # by `corners` on the matrix `mat`. The midpoint of each edge of this square is interpolated # by computing the mean value of the two corners on the edge and the center of the square, and the # displacing it. The displacement is drawn according to `alg.H` and round using `displace` -function _square!(mat, alg::DiamondSquare, round::Int, corners::AbstractVector{Tuple{Int, Int}}) +function _square!(mat, alg::DiamondSquare, round::Int, corners) bottomLeft,bottomRight,topLeft,topRight = corners leftEdge, bottomEdge, topEdge, rightEdge = _edgeMidpointCoordinates(corners) centerPoint = _centerCoordinate(corners) - mat[leftEdge...] = _interpolate(mat, [topLeft,bottomLeft,centerPoint]) + _displace(alg.H, round) - mat[bottomEdge...] = _interpolate(mat, [bottomLeft,bottomRight,centerPoint]) + _displace(alg.H, round) - mat[topEdge...] = _interpolate(mat, [topLeft,topRight,centerPoint]) + _displace(alg.H, round) - mat[rightEdge...] = _interpolate(mat, [topRight,bottomRight,centerPoint]) + _displace(alg.H, round) + mat[leftEdge...] = _interpolate(mat, (topLeft,bottomLeft,centerPoint)) + _displace(alg.H, round) + mat[bottomEdge...] = _interpolate(mat, (bottomLeft,bottomRight,centerPoint)) + _displace(alg.H, round) + mat[topEdge...] = _interpolate(mat, (topLeft,topRight,centerPoint)) + _displace(alg.H, round) + mat[rightEdge...] = _interpolate(mat, (topRight,bottomRight,centerPoint)) + _displace(alg.H, round) end # Runs the square step of the `MidpointDisplacement` algorithm on the square defined # by `corners` on the matrix `mat`. The midpoint of each edge of this square is interpolated # by computing the mean value of the two corners on the edge and the center of the square, and the # displacing it. The displacement is drawn according to `alg.H` and round using `displace` -function _square!(mat, alg::MidpointDisplacement, round::Int, corners::AbstractVector{Tuple{Int, Int}}) +function _square!(mat, alg::MidpointDisplacement, round::Int, corners) bottomLeft,bottomRight,topLeft,topRight = corners leftEdge, bottomEdge, topEdge, rightEdge = _edgeMidpointCoordinates(corners) - mat[leftEdge...] = _interpolate(mat, [topLeft,bottomLeft]) + _displace(alg.H, round) - mat[bottomEdge...] = _interpolate(mat, [bottomLeft,bottomRight]) + _displace(alg.H, round) - mat[topEdge...] = _interpolate(mat, [topLeft,topRight]) + _displace(alg.H, round) - mat[rightEdge...] = _interpolate(mat, [topRight,bottomRight]) + _displace(alg.H, round) + mat[leftEdge...] = _interpolate(mat, (topLeft,bottomLeft)) + _displace(alg.H, round) + mat[bottomEdge...] = _interpolate(mat, (bottomLeft,bottomRight)) + _displace(alg.H, round) + mat[topEdge...] = _interpolate(mat, (topLeft,topRight)) + _displace(alg.H, round) + mat[rightEdge...] = _interpolate(mat, (topRight,bottomRight)) + _displace(alg.H, round) end # Computes the mean of a set of points, represented as a list of indicies to a matrix `mat`. -function _interpolate(mat, points::AbstractVector{Tuple{Int,Int}}) +function _interpolate(mat, points) return mean(mat[pt...] for pt in points) end @@ -169,12 +169,12 @@ end # move from 1.0 to 0 as `round` increases. function _displace(H::Float64, round::Int) σ = 0.5^(round*H) - return(rand(Normal(0, σ))) + return rand(Normal(0, σ)) end # Returns the center coordinate for a square defined by `corners` for the # `DiamondSquare` algorithm. -function _centerCoordinate(corners::AbstractVector{Tuple{Int,Int}}) +function _centerCoordinate(corners) bottomLeft,bottomRight,topLeft,topRight = corners centerX::Int = (_xcoord(bottomLeft)+_xcoord(bottomRight)) ÷ 2 centerY::Int = (_ycoord(topRight)+_ycoord(bottomRight)) ÷ 2 @@ -187,7 +187,7 @@ _xcoord(pt::Tuple{Int,Int}) = pt[1] _ycoord(pt::Tuple{Int,Int}) = pt[2] # Returns an array of midpoints for a square defined by `corners` for the `DiamondSquare` algorithm. -function _edgeMidpointCoordinates(corners::AbstractVector{Tuple{Int,Int}}) +function _edgeMidpointCoordinates(corners) # bottom left, bottom right, top left, top right bottomLeft,bottomRight,topLeft,topRight = corners @@ -196,7 +196,7 @@ function _edgeMidpointCoordinates(corners::AbstractVector{Tuple{Int,Int}}) topEdgeMidpoint::Tuple{Int,Int} = ( (_xcoord(topLeft)+_xcoord(topRight))÷2, _ycoord(topLeft)) rightEdgeMidpoint::Tuple{Int,Int} = ( _xcoord(bottomRight), (_ycoord(bottomRight)+_ycoord(topRight))÷2) - edgeMidpoints = [leftEdgeMidpoint, bottomEdgeMidpoint, topEdgeMidpoint, rightEdgeMidpoint] + edgeMidpoints = (leftEdgeMidpoint, bottomEdgeMidpoint, topEdgeMidpoint, rightEdgeMidpoint) return edgeMidpoints end diff --git a/src/algorithms/distancegradient.jl b/src/algorithms/distancegradient.jl index d2e2b1d..c56b217 100644 --- a/src/algorithms/distancegradient.jl +++ b/src/algorithms/distancegradient.jl @@ -14,11 +14,24 @@ end function _landscape!(mat, alg::DistanceGradient) @assert maximum(alg.sources) <= length(mat) @assert minimum(alg.sources) > 0 - coordinates = Matrix{Float64}(undef, (2, prod(size(mat)))) - for (i, p) in enumerate(Iterators.product(axes(mat)...)) - coordinates[1:2, i] .= p + coordinates = CartesianIndices(mat) + source_coordinates = map(alg.sources) do c + SVector(Tuple(coordinates[c])) + end + idx = Vector{Int}(undef, 1) + dist = Vector{Float64}(undef, 1) + tree = KDTree(source_coordinates) + _write_knn!(mat, dist, idx, tree, coordinates) + return mat +end + +# Function barrier, somehow we lose type stability with this above +function _write_knn!(mat, dist, idx, tree, coordinates) + sortres = false + for i in eachindex(mat) + point = SVector(Tuple(coordinates[i])) + knn_point!(tree, point, sortres, dist, idx, always_false) + mat[i] = dist[1] end - tree = KDTree(coordinates[:,alg.sources]) - mat[:] .= nn(tree, coordinates)[2] return mat end diff --git a/src/algorithms/nncluster.jl b/src/algorithms/nncluster.jl index a52af29..166cec0 100644 --- a/src/algorithms/nncluster.jl +++ b/src/algorithms/nncluster.jl @@ -23,11 +23,21 @@ function _landscape!(mat, alg::NearestNeighborCluster) classify!(mat, [alg.p, 1 - alg.p]) replace!(mat, 2.0 => NaN) clusters, nClusters = label(mat, alg.n) - coordinates = _coordinatematrix(clusters) + coordinates = CartesianIndices(clusters) sources = findall(!isnan, vec(clusters)) - tree = KDTree(coordinates[:,sources]) - clusters[:] .= clusters[sources[nn(tree, coordinates)[1]]] + cluster_coordinates = map(sources) do c + SVector(Tuple(coordinates[c])) + end + idx = Vector{Int}(undef, 1) + dist = Vector{Float64}(undef, 1) + tree = KDTree(cluster_coordinates) randvals = rand(nClusters) - mat .= randvals[Int.(clusters)] - mat + sortres = false + for i in eachindex(mat) + point = SVector(Tuple(coordinates[i])) + knn_point!(tree, point, sortres, dist, idx, always_false) + cluster = clusters[sources[idx[1]]] + mat[i] = randvals[Int(cluster)] + end + return mat end diff --git a/src/algorithms/nnelement.jl b/src/algorithms/nnelement.jl index a5ba8cf..9204481 100644 --- a/src/algorithms/nnelement.jl +++ b/src/algorithms/nnelement.jl @@ -18,26 +18,29 @@ and `k` neighbors. The default is to use three cluster and a single neighbor. end end -function _coordinatematrix(mat) - coordinates = Matrix{Float64}(undef, (2, prod(size(mat)))) - for (i, p) in enumerate(Iterators.product(axes(mat)...)) - coordinates[1:2, i] .= p - end - coordinates -end - function _landscape!(mat, alg::NearestNeighborElement) - fill!(mat, 0.0) clusters = sample(eachindex(mat), alg.n; replace=false) - for (i,n) in enumerate(clusters) - mat[n] = i - end - coordinates = _coordinatematrix(mat) - tree = KDTree(coordinates[:,clusters]) - if alg.k == 1 - mat[:] .= nn(tree, coordinates)[1] + # Preallocate for NearestNeighbors + idx = Vector{Int}(undef, alg.k) + dist = Vector{Float64}(undef, alg.k) + coordinates = CartesianIndices(mat) + cluster_coordinates = map(clusters) do c + SVector(Tuple(coordinates[c])) + end + tree = KDTree(cluster_coordinates) + sortres = false + if alg.k == 1 + for i in eachindex(mat) + point = SVector(Tuple(coordinates[i])) + knn_point!(tree, point, sortres, dist, idx, always_false) + mat[i] = idx[1] + end else - mat[:] .= map(mean, knn(tree, coordinates, alg.k)[1]) + for i in eachindex(mat) + point = SVector(Tuple(coordinates[i])) + knn_point!(tree, point, sortres, dist, idx, always_false) + mat[i] = mean(idx) + end end return mat end diff --git a/src/algorithms/perlinnoise.jl b/src/algorithms/perlinnoise.jl index 3265682..9ce523a 100644 --- a/src/algorithms/perlinnoise.jl +++ b/src/algorithms/perlinnoise.jl @@ -48,10 +48,11 @@ function _landscape!(A, alg::PerlinNoise) # Generate the Perlin noise noise = zeros(eltype(A), (dim, dim)) - meshbuf = Array{eltype(A),3}(undef, dim, dim, 2) + meshbuf1 = Array{eltype(A),2}(undef, dim, dim) + meshbuf2 = Array{eltype(A),2}(undef, dim, dim) nbufs = ntuple(_->Array{eltype(A),2}(undef, dim, dim), 4) for octave in 0:alg.octaves-1 - octave_noise!(noise, meshbuf, nbufs, alg, octave, dim, dim) + octave_noise!(noise, meshbuf1, meshbuf2, nbufs, alg, octave, dim, dim) end # Randomly extract the desired array size @@ -70,7 +71,7 @@ function _landscape!(A, alg::PerlinNoise) end function octave_noise!( - noise, mesh, (n11, n21, n12, n22), alg::PerlinNoise, octave, nrow, ncol + noise, m1, m2, (n11, n21, n12, n22), alg::PerlinNoise, octave, nrow, ncol ) f(t) = @fastmath 6 * t ^ 5 - 15 * t ^ 4 + 10 * t ^ 3 # Wut @@ -78,13 +79,12 @@ function octave_noise!( rp, cp = alg.periods .* alg.lacunarity^(octave) delta = (rp / nrow, cp / ncol) ranges = range(0, rp-delta[1], length=nrow), range(0, cp-delta[2], length=ncol) - _mesh!(mesh, ranges...) - @fastmath mesh .%= 1 + _rem_meshes!(m1, m2, ranges...) # Gradients # This allocates, but the gradients size changes with octave so it needs to # be a very smart in-place `repeat!` or some kind of generator, and the improvement - # may not be than large (~20%). + # may not be that large (~20%). angles = 2pi .* rand(rp + 1, cp + 1) @fastmath gradients = cat(cos.(angles), sin.(angles); dims=3) d = (nrow ÷ rp, ncol ÷ cp) @@ -100,31 +100,28 @@ function octave_noise!( g222 = @view grad[end-nrow+1:end, end-ncol+1:end, 2] # Ramps - m1 = @view mesh[:, :, 1] - m2 = @view mesh[:, :, 2] n11 .= ((m1 .+ m2 ) .* g111 .+ (m1 .+ m2 ) .* g112) n21 .= ((m1 .-1 .+ m2 ) .* g211 .+ (m1 .-1 .+ m2 ) .* g212) n12 .= ((m1 .+ m2 .- 1) .* g121 .+ (m1 .+ m2 .- 1) .* g122) n22 .= ((m1 .-1 .+ m2 .- 1) .* g221 .+ (m1 .-1 .+ m2 .- 1) .* g222) # Interpolation - mesh .= f.(mesh) - t1 = @view mesh[:, :, 1] - t2 = @view mesh[:, :, 2] + m1 .= f.(m1) + m2 .= f.(m2) noise .+= sqrt(2) .* (alg.persistance ^ octave) .* - ((1 .- t2) .* (n11 .* (1 .- t1) .+ t1 .* n21) .+ - t2 .* (n12 .* (1 .- t1) .+ t1 .* n22)) + ((1 .- m2) .* (n11 .* (1 .- m1) .+ m1 .* n21) .+ + m2 .* (n12 .* (1 .- m1) .+ m1 .* n22)) return noise end -function _mesh!(A, x, y) +function _rem_meshes!(m1, m2, x, y) for (i, ival) in enumerate(x), j in 1:length(y) - A[i, j, 1] = ival + @fastmath m1[i, j] = ival % 1 end for i in 1:length(x), (j, jval) in enumerate(y) - A[i, j, 2] = jval + @fastmath m2[i, j] = jval % 1 end - return A + return end function _view_from_square(source, nrow, ncol) diff --git a/src/algorithms/planargradient.jl b/src/algorithms/planargradient.jl index 601b06b..0f1dab2 100644 --- a/src/algorithms/planargradient.jl +++ b/src/algorithms/planargradient.jl @@ -19,6 +19,6 @@ function _landscape!(mat, alg::PlanarGradient) where {IT <: Integer} eastness = sin(deg2rad(alg.direction)) southness = -1cos(deg2rad(alg.direction)) rows, cols = axes(mat) - mat .= southness .* rows .+ (eastness .* cols)' + mat .= collect(rows) .* southness .+ cols' .* eastness return mat end diff --git a/src/algorithms/rectangularcluster.jl b/src/algorithms/rectangularcluster.jl index f2b6ff6..b97eb18 100644 --- a/src/algorithms/rectangularcluster.jl +++ b/src/algorithms/rectangularcluster.jl @@ -19,7 +19,7 @@ end function _landscape!(mat, alg::RectangularCluster) mat .= -1.0 - while minimum(mat) == -1.0 + while any(i -> i === -1.0, mat) width, height = rand(alg.minimum:alg.maximum, 2) row = rand(1:(size(mat,1)-(width-1))) col = rand(1:(size(mat,2)-(height-1))) diff --git a/src/classify.jl b/src/classify.jl index 70d5d14..6d9d0ea 100644 --- a/src/classify.jl +++ b/src/classify.jl @@ -72,9 +72,9 @@ end blend(clusterarray, array, scaling = 1) = blend(clusterarray, [array], [scaling]) const _neighborhoods = Dict( - :rook => ((1, 0), (0, 1)), - :diagonal => ((1, 0), (0, 1), (1, 1)), - :queen => ((1, 0), (0, 1), (1, 1), (1, -1)), + :rook => [(1, 0), (0, 1)], + :diagonal => [(1, 0), (0, 1), (1, 1)], + :queen => [(1, 0), (0, 1), (1, 1), (1, -1)], ) diff --git a/src/landscape.jl b/src/landscape.jl index f1b4a87..f979e8a 100644 --- a/src/landscape.jl +++ b/src/landscape.jl @@ -43,12 +43,11 @@ replaced by `NaN`. function mask!(array::AbstractArray{<:Float64}, maskarray::AbstractArray{<:Bool}) (size(array) == size(maskarray)) || throw(DimensionMismatch("The dimensions of array, $(size(array)), and maskarray, $(size(maskarray)), must match. ")) array[.!maskarray] .= NaN - array + return array end - # Changes the matrix `mat` so that it is between `0` and `1`. function _rescale!(mat) - mat .-= NaNMath.minimum(mat) - mat ./= NaNMath.maximum(mat) + mn, mx = NaNMath.extrema(mat) + mat .= (mat .- mn) ./ (mx - mn) end diff --git a/test/rand.jl b/test/rand.jl index 2b08295..6f8003e 100644 --- a/test/rand.jl +++ b/test/rand.jl @@ -1,5 +1,9 @@ using NeutralLandscapes, Test +@testset "_rescale" begin + @test NeutralLandscapes._rescale!([3.0 -2.0; 0.0 1.0]) == [1.0 0.0; 0.4 0.6] +end + algorithms = ( DiamondSquare(), DiamondSquare(; H=0.2), From 5e1b53600566dbec875f2da9338567df2f5ea106 Mon Sep 17 00:00:00 2001 From: rafaqz Date: Tue, 25 Jan 2022 22:38:31 +0100 Subject: [PATCH 12/31] fix GR docs bug --- docs/make.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/docs/make.jl b/docs/make.jl index 667738c..8288e91 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -1,5 +1,8 @@ using Documenter, NeutralLandscapes +# For GR docs bug +ENV["GKSwstype"] = "100" + makedocs( sitename="Neutral Landscapes", authors="Timothée Poisot", From eda0964dc261afa471990335563282e6223c0a73 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Sat, 2 Apr 2022 10:40:42 -0400 Subject: [PATCH 13/31] :wrench: params/api change --- src/NeutralLandscapes.jl | 7 ++++++- src/updaters/spatial.jl | 8 ++++---- src/updaters/spatiotemporal.jl | 8 ++++---- src/updaters/temporal.jl | 26 +++++++++++++++++++------- src/updaters/update.jl | 17 ++++++++++++++--- 5 files changed, 47 insertions(+), 19 deletions(-) diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index 89139d6..1435c6b 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -44,7 +44,12 @@ include("updaters/update.jl") include("updaters/temporal.jl") include("updaters/spatial.jl") include("updaters/spatiotemporal.jl") -export update, update!, NeutralLandscapeUpdater, TemporallyAutocorrelatedUpdater, SpatiallyAutocorrelatedUpdater, SpatiotemporallyAutocorrelatedUpdater +export update, update! +export NeutralLandscapeUpdater +export TemporallyVariableUpdater +export SpatiallyAutocorrelatedUpdater +export SpatiotemporallyAutocorrelatedUpdater +export rate, variability end # module diff --git a/src/updaters/spatial.jl b/src/updaters/spatial.jl index 695f3e0..0437607 100644 --- a/src/updaters/spatial.jl +++ b/src/updaters/spatial.jl @@ -1,12 +1,12 @@ -@kwdef struct SpatiallyAutocorrelatedUpdater{SU,D,S} <: NeutralLandscapeUpdater +@kwdef struct SpatiallyAutocorrelatedUpdater{SU,R,V} <: NeutralLandscapeUpdater spatialupdater::SU = DiamondSquare(0.5) - direction::D = 0.1 - rate::S = 0.1 + rate::R = 0.1 + variability::V = 0.1 end function _update(sau::SpatiallyAutocorrelatedUpdater, mat) change = rand(spatialupdater(sau), size(mat)) - delta = direction(sau) .+ rate(sau) .* transform(fit(ZScoreTransform, change), change) + delta = rate(sau) .+ variability(sau) .* transform(fit(ZScoreTransform, change), change) mat .+ delta end diff --git a/src/updaters/spatiotemporal.jl b/src/updaters/spatiotemporal.jl index b657301..fa65e17 100644 --- a/src/updaters/spatiotemporal.jl +++ b/src/updaters/spatiotemporal.jl @@ -1,14 +1,14 @@ -@kwdef struct SpatiotemporallyAutocorrelatedUpdater{SU,D,R} <: NeutralLandscapeUpdater +@kwdef struct SpatiotemporallyAutocorrelatedUpdater{SU,R,V} <: NeutralLandscapeUpdater spatialupdater::SU = DiamondSquare(0.1) - direction::D = 0.1 rate::R = 0.1 + variability::V = 0.1 end function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) change = rand(spatialupdater(stau), size(mat)) - temporalshift = rand(Normal(), size(mat)) + temporalshift = rand(Normal(0, variability(stau)), size(mat)) z = transform(fit(ZScoreTransform, change), change) - delta = direction(stau) .+ (rate(stau) .* z) .+ (rate(stau) .* temporalshift) + delta = rate(stau) .+ z .+ temporalshift mat .+ delta end \ No newline at end of file diff --git a/src/updaters/temporal.jl b/src/updaters/temporal.jl index 89a3861..d2c6533 100644 --- a/src/updaters/temporal.jl +++ b/src/updaters/temporal.jl @@ -1,11 +1,23 @@ +""" +TemporallyVariableUpdater{D,S} <: NeutralLandscapeUpdater -@kwdef struct TemporallyAutocorrelatedUpdater{D,S} <: NeutralLandscapeUpdater - spatialupdater = missing - direction::D = 0.1 - rate::S = 0.1 +A `NeutralLandscapeUpdater` that has a prescribed + +""" +@kwdef struct TemporallyVariableUpdater{D,R,V} <: NeutralLandscapeUpdater + spatialupdater::D = missing + rate::R = 0.1 + variability::V = 0.1 end -function _update(sau::TemporallyAutocorrelatedUpdater, mat) - U = rand(Normal(), size(mat)) - Δ = direction(sau) .+ rate(sau) .* U + +""" +_update(tvu::TemporallyVariableUpdater, mat) + +Updates `mat` using temporally autocorrelated change, +using the direction and rate parameters from `tvu`. +""" +function _update(tvu::TemporallyVariableUpdater, mat) + U = rand(Normal(0, variability(tvu)), size(mat)) + Δ = rate(tvu) .+ U mat .+ Δ end diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 0e9f4ca..eed63d6 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -15,8 +15,19 @@ abstract type NeutralLandscapeUpdater end which defines the expected change in any cell per timestep. """ rate(up::NeutralLandscapeUpdater) = up.rate + +""" + All `NeutralLandscapeUpdater`'s have a `spatialupdater` field + which is either a `NeutralLandscapeMaker`, or `Missing` (in the case + of temporally correlated updaters). +""" spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater -direction(up::NeutralLandscapeUpdater) = up.direction + +""" + TODO rate and direction should have different names to clarify what they actually are + +""" +variability(up::NeutralLandscapeUpdater) = up.variability function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} @@ -24,13 +35,13 @@ function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} end function update(updater::T, mat, n::I) where {T<:NeutralLandscapeUpdater, I<:Integer} sequence = [zeros(size(mat)) for _ in 1:n] - sequence[begin] = mat + sequence[begin] .= mat for i in 2:n sequence[i] = _update(updater, sequence[i-1]) end sequence end function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} - mat[:] = _update(updater, mat) + mat .= _update(updater, mat) end From 6510d39f3529e1ab2bbf5391f37346b4c6a9876e Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 5 Apr 2022 11:44:44 -0400 Subject: [PATCH 14/31] normalization of time-series of landscapes --- src/NeutralLandscapes.jl | 3 +++ src/updaters/update.jl | 20 ++++++++++++++++++++ 2 files changed, 23 insertions(+) diff --git a/src/NeutralLandscapes.jl b/src/NeutralLandscapes.jl index 880db0a..0aad384 100644 --- a/src/NeutralLandscapes.jl +++ b/src/NeutralLandscapes.jl @@ -50,7 +50,10 @@ export TemporallyVariableUpdater export SpatiallyAutocorrelatedUpdater export SpatiotemporallyAutocorrelatedUpdater export rate, variability +export normalize end # module + + diff --git a/src/updaters/update.jl b/src/updaters/update.jl index eed63d6..0ef42ce 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -29,6 +29,25 @@ spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater """ variability(up::NeutralLandscapeUpdater) = up.variability +""" + normalize +""" +function normalize(mats::Vector{M}) where {M<:AbstractMatrix} + + # We want to normalize such that all values fall on 0 to 1 scale, + # where the min value seen across all matrices in `vec` + + mins, maxs = findmin.(mats), findmax.(mats) + totalmin, totalmax = findmin([x[1] for x in mins])[1], findmax([x[1] for x in maxs])[1] + + returnmats = copy(mats) + + for (i,mat) in enumerate(mats) + returnmats[i] = (mat .+ totalmin) ./ totalmax + end + + return returnmats +end function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) @@ -45,3 +64,4 @@ function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} mat .= _update(updater, mat) end + From b35e036af7623ac2ec554c5737a2672c7cab310f Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 5 Apr 2022 11:53:44 -0400 Subject: [PATCH 15/31] :bug: fixed bug in normalize --- src/updaters/update.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 0ef42ce..2e4a1aa 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -39,11 +39,11 @@ function normalize(mats::Vector{M}) where {M<:AbstractMatrix} mins, maxs = findmin.(mats), findmax.(mats) totalmin, totalmax = findmin([x[1] for x in mins])[1], findmax([x[1] for x in maxs])[1] - + @show totalmin, totalmax returnmats = copy(mats) for (i,mat) in enumerate(mats) - returnmats[i] = (mat .+ totalmin) ./ totalmax + returnmats[i] = (mat .- totalmin) ./ (totalmax - totalmin) end return returnmats From 6ef5bba99aa5b0f8378426fc761917c648c597e6 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 5 Apr 2022 12:02:35 -0400 Subject: [PATCH 16/31] :bug: fixed bug in stau implementation --- src/updaters/spatiotemporal.jl | 2 +- src/updaters/update.jl | 29 ++++++++++++++++++----------- 2 files changed, 19 insertions(+), 12 deletions(-) diff --git a/src/updaters/spatiotemporal.jl b/src/updaters/spatiotemporal.jl index fa65e17..79eefa5 100644 --- a/src/updaters/spatiotemporal.jl +++ b/src/updaters/spatiotemporal.jl @@ -9,6 +9,6 @@ function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) change = rand(spatialupdater(stau), size(mat)) temporalshift = rand(Normal(0, variability(stau)), size(mat)) z = transform(fit(ZScoreTransform, change), change) - delta = rate(stau) .+ z .+ temporalshift + delta = rate(stau) .+ variability(stau) * z .+ temporalshift mat .+ delta end \ No newline at end of file diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 2e4a1aa..f9dac78 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -12,7 +12,7 @@ abstract type NeutralLandscapeUpdater end """ All `NeutralLandscapeUpdater`s have a field `rate` - which defines the expected change in any cell per timestep. + which defines the expected (or mean) change across all cells per timestep. """ rate(up::NeutralLandscapeUpdater) = up.rate @@ -24,34 +24,40 @@ rate(up::NeutralLandscapeUpdater) = up.rate spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater """ - TODO rate and direction should have different names to clarify what they actually are - + variability(up::NeutralLandscapeUpdater) + + Returns the `variability` of a given `NeutralLandscapeUpdater`. + The variability of an updater is how much temporal variation there + will be in a generated time-series of landscapes. """ variability(up::NeutralLandscapeUpdater) = up.variability """ - normalize -""" -function normalize(mats::Vector{M}) where {M<:AbstractMatrix} + normalize(mats::Vector{M}) - # We want to normalize such that all values fall on 0 to 1 scale, - # where the min value seen across all matrices in `vec` + Normalizes a vector of neutral landscapes `mats` such that all + values between 0 and 1. Note that this does not preserve the + `rate` parameter for a given `NeutralLandscapeUpdater`, and instead + rescales it proportional to the difference between the total maximum + and total minimum across all `mats`. +""" +function normalize(mats::Vector{M}) where {M<:AbstractMatrix} mins, maxs = findmin.(mats), findmax.(mats) totalmin, totalmax = findmin([x[1] for x in mins])[1], findmax([x[1] for x in maxs])[1] - @show totalmin, totalmax returnmats = copy(mats) - for (i,mat) in enumerate(mats) returnmats[i] = (mat .- totalmin) ./ (totalmax - totalmin) end - return returnmats end + + function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) end + function update(updater::T, mat, n::I) where {T<:NeutralLandscapeUpdater, I<:Integer} sequence = [zeros(size(mat)) for _ in 1:n] sequence[begin] .= mat @@ -60,6 +66,7 @@ function update(updater::T, mat, n::I) where {T<:NeutralLandscapeUpdater, I<:Int end sequence end + function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} mat .= _update(updater, mat) end From 2605e02d4fe3e5f4d4b9118dc4522df9128fb413 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 5 Apr 2022 12:14:10 -0400 Subject: [PATCH 17/31] :pencil: all docstrings exist --- src/updaters/spatial.jl | 19 +++++++++++++++++++ src/updaters/spatiotemporal.jl | 22 ++++++++++++++++++++++ src/updaters/temporal.jl | 5 +++-- src/updaters/update.jl | 18 ++++++++++++++++++ 4 files changed, 62 insertions(+), 2 deletions(-) diff --git a/src/updaters/spatial.jl b/src/updaters/spatial.jl index 0437607..21e50bb 100644 --- a/src/updaters/spatial.jl +++ b/src/updaters/spatial.jl @@ -1,10 +1,29 @@ +""" +SpatiallyAutocorrelatedUpdater{SU,R,V} +A `NeutralLandscapeUpdater` that has a prescribed level of +spatial variation (`variability`) and rate of change (`rate`), +and where the spatial distribution of this change is proportional +to a neutral landscape generated with `spatialupdater` at every time +step. + +TODO: make it possible to fix a given spatial updater at each timestep. +""" @kwdef struct SpatiallyAutocorrelatedUpdater{SU,R,V} <: NeutralLandscapeUpdater spatialupdater::SU = DiamondSquare(0.5) rate::R = 0.1 variability::V = 0.1 end + +""" +_update(sau::SpatiallyAutocorrelatedUpdater, mat) + +Updates `mat` using spatially autocorrelated change, +using the direction, rate, and spatialupdater parameters from `sau" + +TODO: doesn't necessarily have to be a ZScoreTransform, could be arbitrary argument +""" function _update(sau::SpatiallyAutocorrelatedUpdater, mat) change = rand(spatialupdater(sau), size(mat)) delta = rate(sau) .+ variability(sau) .* transform(fit(ZScoreTransform, change), change) diff --git a/src/updaters/spatiotemporal.jl b/src/updaters/spatiotemporal.jl index 79eefa5..5705363 100644 --- a/src/updaters/spatiotemporal.jl +++ b/src/updaters/spatiotemporal.jl @@ -1,10 +1,32 @@ +""" +SpatiotemporallyAutocorrelatedUpdater{SU,R,V} + +A `NeutralLandscapeUpdater` that has a prescribed level of +spatial and temporal variation (`variability`) and rate of change (`rate`), +and where the spatial distribution of this change is proportional +to a neutral landscape generated with `spatialupdater` at every time +step. + +# TODO perhaps spatial and temporal should each have their own variability param +""" @kwdef struct SpatiotemporallyAutocorrelatedUpdater{SU,R,V} <: NeutralLandscapeUpdater spatialupdater::SU = DiamondSquare(0.1) rate::R = 0.1 variability::V = 0.1 end + + +""" +_update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) + +Updates `mat` using temporally autocorrelated change, +using the direction, rate, and spatialupdater parameters from `stau" + +TODO: doesn't necessarily have to be a Normal distribution or ZScoreTransform, +could be arbitrary argument +""" function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) change = rand(spatialupdater(stau), size(mat)) temporalshift = rand(Normal(0, variability(stau)), size(mat)) diff --git a/src/updaters/temporal.jl b/src/updaters/temporal.jl index d2c6533..c3b938f 100644 --- a/src/updaters/temporal.jl +++ b/src/updaters/temporal.jl @@ -1,8 +1,9 @@ """ TemporallyVariableUpdater{D,S} <: NeutralLandscapeUpdater -A `NeutralLandscapeUpdater` that has a prescribed - +A `NeutralLandscapeUpdater` that has a prescribed level of +temporal variation (`variability`) and rate of change (`rate`), +but no spatial correlation in where change is distributed. """ @kwdef struct TemporallyVariableUpdater{D,R,V} <: NeutralLandscapeUpdater spatialupdater::D = missing diff --git a/src/updaters/update.jl b/src/updaters/update.jl index f9dac78..0cca000 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -53,11 +53,23 @@ function normalize(mats::Vector{M}) where {M<:AbstractMatrix} end +""" + update(updater::T, mat) + Returns one-timestep applied to `mat` based on the `NeutralLandscapeUpdater` + provided (`updater`). +""" function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) end + +""" +update(updater::T, mat, n::I) + +Returns a sequence of length `n` where the original neutral landscape +`mat` is updated by the `NeutralLandscapeUpdater` `update` for `n` timesteps. +""" function update(updater::T, mat, n::I) where {T<:NeutralLandscapeUpdater, I<:Integer} sequence = [zeros(size(mat)) for _ in 1:n] sequence[begin] .= mat @@ -67,6 +79,12 @@ function update(updater::T, mat, n::I) where {T<:NeutralLandscapeUpdater, I<:Int sequence end + +""" +update!(updater::T, mat) + +Updates a landscape `mat` in-place by directly mutating `mat`. +""" function update!(updater::T, mat) where {T<:NeutralLandscapeUpdater} mat .= _update(updater, mat) end From 0cd8b69bbbab318ac956e2f75cb4e027dda37c77 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 5 Apr 2022 12:20:17 -0400 Subject: [PATCH 18/31] added test file for updaters --- test/runtests.jl | 1 + test/updaters.jl | 7 +++++++ 2 files changed, 8 insertions(+) create mode 100644 test/updaters.jl diff --git a/test/runtests.jl b/test/runtests.jl index 98d5327..dcabf80 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,4 +1,5 @@ using Test, SafeTestsets +@time @safetestset "updaters" begin include("updaters.jl") end @time @safetestset "planar gradient" begin include("planar.jl") end @time @safetestset "rand" begin include("rand.jl") end diff --git a/test/updaters.jl b/test/updaters.jl new file mode 100644 index 0000000..0e2be47 --- /dev/null +++ b/test/updaters.jl @@ -0,0 +1,7 @@ +using NeutralLandscapes +using Test + +# test it's running +@test TemporallyVariableUpdater() != π +@test SpatiallyAutocorrelatedUpdater() != π +@test SpatiotemporallyAutocorrelatedUpdater() != π \ No newline at end of file From 0256a1b2b0663ea01d52cbbbb3b77141b474c6aa Mon Sep 17 00:00:00 2001 From: michael catchen Date: Wed, 6 Apr 2022 09:00:07 -0400 Subject: [PATCH 19/31] more tests --- src/updaters/temporal.jl | 6 ++++++ test/updaters.jl | 36 +++++++++++++++++++++++++++++++++++- 2 files changed, 41 insertions(+), 1 deletion(-) diff --git a/src/updaters/temporal.jl b/src/updaters/temporal.jl index c3b938f..0ef47c3 100644 --- a/src/updaters/temporal.jl +++ b/src/updaters/temporal.jl @@ -16,6 +16,12 @@ _update(tvu::TemporallyVariableUpdater, mat) Updates `mat` using temporally autocorrelated change, using the direction and rate parameters from `tvu`. + + +#TODO this doesn't have to be a Normal distribution, +#could be arbitrary distribution that is continuous +#and can have mean 0 (or that can be transformed to +# have mean 0) """ function _update(tvu::TemporallyVariableUpdater, mat) U = rand(Normal(0, variability(tvu)), size(mat)) diff --git a/test/updaters.jl b/test/updaters.jl index 0e2be47..7f51733 100644 --- a/test/updaters.jl +++ b/test/updaters.jl @@ -4,4 +4,38 @@ using Test # test it's running @test TemporallyVariableUpdater() != π @test SpatiallyAutocorrelatedUpdater() != π -@test SpatiotemporallyAutocorrelatedUpdater() != π \ No newline at end of file +@test SpatiotemporallyAutocorrelatedUpdater() != π + +function testupdaters(model) + updater = model() + + # Test defaults + @test rate(updater) == 0.1 + @test variability(updater) == 0.1 + + # Test kwargs + updater = model(rate = 1.0, variability=0.05, spatialupdater=MidpointDisplacement(0.5)) + @test rate(updater) == 1.0 + @test variability(updater) == 0.05 + @test typeof(updater.spatialupdater) <: NeutralLandscapeMaker + @test updater.spatialupdater == MidpointDisplacement(0.5) + + # Test updating + env = rand(MidpointDisplacement(0.5), 50, 50) + newenv = update(updater, env) + @test env != newenv + + oldenv = deepcopy(env) + update!(updater, env) + @test env != oldenv +end + + + +models = [ + TemporallyVariableUpdater, + SpatiallyAutocorrelatedUpdater, + SpatiotemporallyAutocorrelatedUpdater +] + +testupdaters.(models) From 0710b754bb2c2b4f899b987b3f7736bff92e116c Mon Sep 17 00:00:00 2001 From: michael catchen Date: Mon, 11 Apr 2022 13:43:58 -0400 Subject: [PATCH 20/31] fix doc indents --- src/updaters/update.jl | 38 +++++++++++++++++++------------------- 1 file changed, 19 insertions(+), 19 deletions(-) diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 0cca000..3067b0b 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -1,25 +1,25 @@ """ - This file contains methods for simulating change in a - landscape with expected amounts of spatial and temporal - autocorrelation. +This file contains methods for simulating change in a +landscape with expected amounts of spatial and temporal +autocorrelation. """ """ - NeutralLandscapeUpdater is an abstract type for methods - for updating a landscape matrix +NeutralLandscapeUpdater is an abstract type for methods +for updating a landscape matrix """ abstract type NeutralLandscapeUpdater end """ - All `NeutralLandscapeUpdater`s have a field `rate` - which defines the expected (or mean) change across all cells per timestep. +All `NeutralLandscapeUpdater`s have a field `rate` +which defines the expected (or mean) change across all cells per timestep. """ rate(up::NeutralLandscapeUpdater) = up.rate """ - All `NeutralLandscapeUpdater`'s have a `spatialupdater` field - which is either a `NeutralLandscapeMaker`, or `Missing` (in the case - of temporally correlated updaters). +All `NeutralLandscapeUpdater`'s have a `spatialupdater` field +which is either a `NeutralLandscapeMaker`, or `Missing` (in the case +of temporally correlated updaters). """ spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater @@ -33,14 +33,14 @@ spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater variability(up::NeutralLandscapeUpdater) = up.variability """ - normalize(mats::Vector{M}) +normalize(mats::Vector{M}) - Normalizes a vector of neutral landscapes `mats` such that all - values between 0 and 1. Note that this does not preserve the - `rate` parameter for a given `NeutralLandscapeUpdater`, and instead - rescales it proportional to the difference between the total maximum - and total minimum across all `mats`. +Normalizes a vector of neutral landscapes `mats` such that all +values between 0 and 1. Note that this does not preserve the +`rate` parameter for a given `NeutralLandscapeUpdater`, and instead +rescales it proportional to the difference between the total maximum +and total minimum across all `mats`. """ function normalize(mats::Vector{M}) where {M<:AbstractMatrix} mins, maxs = findmin.(mats), findmax.(mats) @@ -54,10 +54,10 @@ end """ - update(updater::T, mat) +update(updater::T, mat) - Returns one-timestep applied to `mat` based on the `NeutralLandscapeUpdater` - provided (`updater`). +Returns one-timestep applied to `mat` based on the `NeutralLandscapeUpdater` +provided (`updater`). """ function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} _update(updater, mat) From 9523413a063ea3fc075f46cc8c2483c4908dbd25 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Mon, 11 Apr 2022 16:47:01 -0400 Subject: [PATCH 21/31] docstring fix --- src/updaters/update.jl | 6 ------ 1 file changed, 6 deletions(-) diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 3067b0b..e09d4c0 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -1,9 +1,3 @@ -""" -This file contains methods for simulating change in a -landscape with expected amounts of spatial and temporal -autocorrelation. -""" - """ NeutralLandscapeUpdater is an abstract type for methods for updating a landscape matrix From ca61f735a3b1de870150b26825fa0e001302c574 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Mon, 11 Apr 2022 16:53:57 -0400 Subject: [PATCH 22/31] more docstinrg fixes --- src/updaters/update.jl | 16 +++++++++++----- 1 file changed, 11 insertions(+), 5 deletions(-) diff --git a/src/updaters/update.jl b/src/updaters/update.jl index e09d4c0..1ee441e 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -1,16 +1,22 @@ """ +NeutralLandscapeUpdater + NeutralLandscapeUpdater is an abstract type for methods for updating a landscape matrix """ abstract type NeutralLandscapeUpdater end """ +spatialupdater(up::NeutralLandscapeUpdater) + All `NeutralLandscapeUpdater`s have a field `rate` which defines the expected (or mean) change across all cells per timestep. """ rate(up::NeutralLandscapeUpdater) = up.rate """ +spatialupdater(up::NeutralLandscapeUpdater) + All `NeutralLandscapeUpdater`'s have a `spatialupdater` field which is either a `NeutralLandscapeMaker`, or `Missing` (in the case of temporally correlated updaters). @@ -18,11 +24,11 @@ of temporally correlated updaters). spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater """ - variability(up::NeutralLandscapeUpdater) - - Returns the `variability` of a given `NeutralLandscapeUpdater`. - The variability of an updater is how much temporal variation there - will be in a generated time-series of landscapes. +variability(up::NeutralLandscapeUpdater) + +Returns the `variability` of a given `NeutralLandscapeUpdater`. +The variability of an updater is how much temporal variation there +will be in a generated time-series of landscapes. """ variability(up::NeutralLandscapeUpdater) = up.variability From cc8a4851187dd9492eaf8cb3c3b4c4953a10aec9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Timoth=C3=A9e=20Poisot?= Date: Mon, 11 Apr 2022 20:41:50 -0400 Subject: [PATCH 23/31] Bump NaNMath to 1.0, closes #60 --- Project.toml | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index f79b0ed..adfb6f6 100644 --- a/Project.toml +++ b/Project.toml @@ -8,7 +8,6 @@ DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" NearestNeighbors = "b8a86587-4115-5ab1-83bc-aa920d37bbce" -ProfileView = "c46f51b8-102a-5cf2-8d2c-8597cb0e0da7" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" @@ -16,7 +15,7 @@ StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" [compat] DataStructures = "0.18" Distributions = "0.24, 0.25" -NaNMath = "0.3" +NaNMath = "1.0" NearestNeighbors = "0.4" StatsBase = "0.33" julia = "1.4" From 869068eb75dd35b9b8bce873892ebe3cb5377329 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Timoth=C3=A9e=20Poisot?= Date: Mon, 11 Apr 2022 20:42:04 -0400 Subject: [PATCH 24/31] fix the indentation in docstrings --- src/updaters/spatial.jl | 11 +++---- src/updaters/spatiotemporal.jl | 13 ++++----- src/updaters/temporal.jl | 22 +++++++------- src/updaters/update.jl | 52 ++++++++++++++++------------------ 4 files changed, 47 insertions(+), 51 deletions(-) diff --git a/src/updaters/spatial.jl b/src/updaters/spatial.jl index 21e50bb..523bb74 100644 --- a/src/updaters/spatial.jl +++ b/src/updaters/spatial.jl @@ -1,6 +1,6 @@ """ -SpatiallyAutocorrelatedUpdater{SU,R,V} + SpatiallyAutocorrelatedUpdater{SU,R,V} A `NeutralLandscapeUpdater` that has a prescribed level of spatial variation (`variability`) and rate of change (`rate`), @@ -17,12 +17,13 @@ TODO: make it possible to fix a given spatial updater at each timestep. end """ -_update(sau::SpatiallyAutocorrelatedUpdater, mat) + _update(sau::SpatiallyAutocorrelatedUpdater, mat; transform=ZScoreTransform) -Updates `mat` using spatially autocorrelated change, -using the direction, rate, and spatialupdater parameters from `sau" +Updates `mat` using spatially autocorrelated change, using the direction, rate, +and spatial updater parameters from `sau`. -TODO: doesn't necessarily have to be a ZScoreTransform, could be arbitrary argument +TODO: doesn't necessarily have to be a ZScoreTransform, could be arbitrary +argument """ function _update(sau::SpatiallyAutocorrelatedUpdater, mat) change = rand(spatialupdater(sau), size(mat)) diff --git a/src/updaters/spatiotemporal.jl b/src/updaters/spatiotemporal.jl index 5705363..3f115d9 100644 --- a/src/updaters/spatiotemporal.jl +++ b/src/updaters/spatiotemporal.jl @@ -1,7 +1,6 @@ """ -SpatiotemporallyAutocorrelatedUpdater{SU,R,V} - + SpatiotemporallyAutocorrelatedUpdater{SU,R,V} A `NeutralLandscapeUpdater` that has a prescribed level of spatial and temporal variation (`variability`) and rate of change (`rate`), @@ -9,7 +8,7 @@ and where the spatial distribution of this change is proportional to a neutral landscape generated with `spatialupdater` at every time step. -# TODO perhaps spatial and temporal should each have their own variability param +TODO: perhaps spatial and temporal should each have their own variability param """ @kwdef struct SpatiotemporallyAutocorrelatedUpdater{SU,R,V} <: NeutralLandscapeUpdater spatialupdater::SU = DiamondSquare(0.1) @@ -19,12 +18,12 @@ end """ -_update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) + _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) -Updates `mat` using temporally autocorrelated change, -using the direction, rate, and spatialupdater parameters from `stau" +Updates `mat` using temporally autocorrelated change, using the direction, rate, +and spatialupdater parameters from `stau`. -TODO: doesn't necessarily have to be a Normal distribution or ZScoreTransform, +TODO: doesn't necessarily have to be a Normal distribution or ZScoreTransform, could be arbitrary argument """ function _update(stau::SpatiotemporallyAutocorrelatedUpdater, mat) diff --git a/src/updaters/temporal.jl b/src/updaters/temporal.jl index 0ef47c3..6f4f8f8 100644 --- a/src/updaters/temporal.jl +++ b/src/updaters/temporal.jl @@ -1,9 +1,9 @@ """ -TemporallyVariableUpdater{D,S} <: NeutralLandscapeUpdater + TemporallyVariableUpdater{D,S} <: NeutralLandscapeUpdater -A `NeutralLandscapeUpdater` that has a prescribed level of -temporal variation (`variability`) and rate of change (`rate`), -but no spatial correlation in where change is distributed. +A `NeutralLandscapeUpdater` that has a prescribed level of temporal variation +(`variability`) and rate of change (`rate`), but no spatial correlation in where +change is distributed. """ @kwdef struct TemporallyVariableUpdater{D,R,V} <: NeutralLandscapeUpdater spatialupdater::D = missing @@ -12,16 +12,14 @@ but no spatial correlation in where change is distributed. end """ -_update(tvu::TemporallyVariableUpdater, mat) + _update(tvu::TemporallyVariableUpdater, mat) -Updates `mat` using temporally autocorrelated change, -using the direction and rate parameters from `tvu`. +Updates `mat` using temporally autocorrelated change, using the direction and +rate parameters from `tvu`. - -#TODO this doesn't have to be a Normal distribution, -#could be arbitrary distribution that is continuous -#and can have mean 0 (or that can be transformed to -# have mean 0) +TODO: this doesn't have to be a Normal distribution, could be arbitrary +distribution that is continuous and can have mean 0 (or that can be transformed +to have mean 0) """ function _update(tvu::TemporallyVariableUpdater, mat) U = rand(Normal(0, variability(tvu)), size(mat)) diff --git a/src/updaters/update.jl b/src/updaters/update.jl index 1ee441e..0434d92 100644 --- a/src/updaters/update.jl +++ b/src/updaters/update.jl @@ -1,46 +1,44 @@ """ -NeutralLandscapeUpdater + NeutralLandscapeUpdater -NeutralLandscapeUpdater is an abstract type for methods -for updating a landscape matrix +NeutralLandscapeUpdater is an abstract type for methods for updating a landscape +matrix """ abstract type NeutralLandscapeUpdater end """ -spatialupdater(up::NeutralLandscapeUpdater) + spatialupdater(up::NeutralLandscapeUpdater) -All `NeutralLandscapeUpdater`s have a field `rate` -which defines the expected (or mean) change across all cells per timestep. +All `NeutralLandscapeUpdater`s have a field `rate` which defines the expected +(or mean) change across all cells per timestep. """ rate(up::NeutralLandscapeUpdater) = up.rate """ -spatialupdater(up::NeutralLandscapeUpdater) + spatialupdater(up::NeutralLandscapeUpdater) -All `NeutralLandscapeUpdater`'s have a `spatialupdater` field -which is either a `NeutralLandscapeMaker`, or `Missing` (in the case -of temporally correlated updaters). +All `NeutralLandscapeUpdater`'s have a `spatialupdater` field which is either a +`NeutralLandscapeMaker`, or `Missing` (in the case of temporally correlated +updaters). """ spatialupdater(up::NeutralLandscapeUpdater) = up.spatialupdater """ -variability(up::NeutralLandscapeUpdater) + variability(up::NeutralLandscapeUpdater) -Returns the `variability` of a given `NeutralLandscapeUpdater`. -The variability of an updater is how much temporal variation there -will be in a generated time-series of landscapes. +Returns the `variability` of a given `NeutralLandscapeUpdater`. The variability +of an updater is how much temporal variation there will be in a generated +time-series of landscapes. """ variability(up::NeutralLandscapeUpdater) = up.variability """ -normalize(mats::Vector{M}) + normalize(mats::Vector{M}) - -Normalizes a vector of neutral landscapes `mats` such that all -values between 0 and 1. Note that this does not preserve the -`rate` parameter for a given `NeutralLandscapeUpdater`, and instead -rescales it proportional to the difference between the total maximum -and total minimum across all `mats`. +Normalizes a vector of neutral landscapes `mats` such that all values between 0 +and 1. Note that this does not preserve the `rate` parameter for a given +`NeutralLandscapeUpdater`, and instead rescales it proportional to the +difference between the total maximum and total minimum across all `mats`. """ function normalize(mats::Vector{M}) where {M<:AbstractMatrix} mins, maxs = findmin.(mats), findmax.(mats) @@ -54,9 +52,9 @@ end """ -update(updater::T, mat) + update(updater::T, mat) -Returns one-timestep applied to `mat` based on the `NeutralLandscapeUpdater` +Returns one-timestep applied to `mat` based on the `NeutralLandscapeUpdater` provided (`updater`). """ function update(updater::T, mat) where {T<:NeutralLandscapeUpdater} @@ -65,10 +63,10 @@ end """ -update(updater::T, mat, n::I) + update(updater::T, mat, n::I) -Returns a sequence of length `n` where the original neutral landscape -`mat` is updated by the `NeutralLandscapeUpdater` `update` for `n` timesteps. +Returns a sequence of length `n` where the original neutral landscape `mat` is +updated by the `NeutralLandscapeUpdater` `update` for `n` timesteps. """ function update(updater::T, mat, n::I) where {T<:NeutralLandscapeUpdater, I<:Integer} sequence = [zeros(size(mat)) for _ in 1:n] @@ -81,7 +79,7 @@ end """ -update!(updater::T, mat) + update!(updater::T, mat) Updates a landscape `mat` in-place by directly mutating `mat`. """ From 9dc7edcf9f080589dd1e68959f20938e8e17b963 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Timoth=C3=A9e=20Poisot?= Date: Mon, 11 Apr 2022 20:45:15 -0400 Subject: [PATCH 25/31] update NaNMath possible versions --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index adfb6f6..de495b1 100644 --- a/Project.toml +++ b/Project.toml @@ -15,7 +15,7 @@ StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" [compat] DataStructures = "0.18" Distributions = "0.24, 0.25" -NaNMath = "1.0" +NaNMath = "0.3, 1.0" NearestNeighbors = "0.4" StatsBase = "0.33" julia = "1.4" From 7e6a769c948e5c5de18880423a345a5b6187a47b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Timoth=C3=A9e=20Poisot?= Date: Mon, 11 Apr 2022 20:45:45 -0400 Subject: [PATCH 26/31] add Julia v1.6 and v1.7 to CI --- .github/workflows/CI.yml | 2 ++ 1 file changed, 2 insertions(+) diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml index d6eaf94..3851b14 100644 --- a/.github/workflows/CI.yml +++ b/.github/workflows/CI.yml @@ -18,6 +18,8 @@ jobs: version: - '1.4' - '1.5' + - '1.6' + - '1.7' os: - ubuntu-latest - macOS-latest From 811dd081d49cc47a07cddf6cdf648fe6724242bb Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 12 Apr 2022 11:24:42 -0400 Subject: [PATCH 27/31] docs update --- docs/src/index.md | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/docs/src/index.md b/docs/src/index.md index 4baa78e..fbda44e 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -27,6 +27,22 @@ rand rand! ``` +## Temporal Change +```@docs + NeutralLandscapeUpdater + TemporallyVariableUpdater + SpatiallyAutocorrelatedUpdater + SpatiotemporallyAutocorrelatedUpdater + update + update! + spatialupdater, + variability, + rate + _update + normalize +``` + + ## Other functions ```@docs From 4f636cfd0aa14f82412b25674a5f58eee9bee3c6 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 12 Apr 2022 11:34:01 -0400 Subject: [PATCH 28/31] internal _updater remove --- docs/src/index.md | 1 - 1 file changed, 1 deletion(-) diff --git a/docs/src/index.md b/docs/src/index.md index fbda44e..db20981 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -38,7 +38,6 @@ rand! spatialupdater, variability, rate - _update normalize ``` From 193ae0f9d39534d08b751f34e471976bd02425ec Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 12 Apr 2022 11:46:33 -0400 Subject: [PATCH 29/31] remove internal function from docs --- docs/src/index.md | 3 --- 1 file changed, 3 deletions(-) diff --git a/docs/src/index.md b/docs/src/index.md index db20981..60f54e2 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -35,9 +35,6 @@ rand! SpatiotemporallyAutocorrelatedUpdater update update! - spatialupdater, - variability, - rate normalize ``` From 539aaa03771f2b9569fab274d10b34044602648e Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 12 Apr 2022 12:01:46 -0400 Subject: [PATCH 30/31] trying to get documenter to recognize internal funcs --- docs/src/index.md | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/docs/src/index.md b/docs/src/index.md index 60f54e2..a16307f 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -36,6 +36,10 @@ rand! update update! normalize + NeutralLandscapes.rate + NeutralLandscapes.variability + NeutralLandscapes.spatialupdater + NeutralLandscapes._update ``` From 7058ba3330341fc353eeaf1441d0291ea32597e5 Mon Sep 17 00:00:00 2001 From: michael catchen Date: Tue, 12 Apr 2022 12:34:17 -0400 Subject: [PATCH 31/31] bump version --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index de495b1..97aedec 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "NeutralLandscapes" uuid = "71847384-8354-4223-ac08-659a5128069f" authors = ["Timothée Poisot ", "Michael Krabbe Borregaard ", "Michael David Catchen ", "Rafael Schouten ", "Virgile Baudrot "] -version = "0.0.4" +version = "0.1.0" [deps] DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"