Using a custom dispersal density functon for indviduals

[Sherif Negm]

 I have a custom density function for dispersal in a continuous habitat. How can i write it down, so individuals can disperse using it. According to the manual, it only showed how individuals can disperse using preset functions like normal distribution density. 

[Peter Ralph]

Hello, Sherif! Basically, you need to figure out how to draw a random sample from your density function using the tools we have available. The most likely options are (a) sample using an implemented kernel and rejection sample; (b) draw a random direction and a random distance by putting a Uniform(0,1) through the inverse CDF of the distance for your kernel - but, it's a big topic.

--peter

---

From: slim-d...@googlegroups.com <slim-d...@googlegroups.com> on behalf of Sherif Negm <sheref...@gmail.com>
Sent: Tuesday, April 16, 2024 8:02 AM
To: slim-discuss <slim-d...@googlegroups.com>
Subject: Using a custom dispersal density functon for indviduals

 

 I have a custom density function for dispersal in a continuous habitat. How can i write it down, so individuals can disperse using it. According to the manual, it only showed how individuals can disperse using preset functions like normal distribution density. 

--
SLiM forward genetic simulation: http://messerlab.org/slim/
---
You received this message because you are subscribed to the Google Groups "slim-discuss" group.
To unsubscribe from this group and stop receiving emails from it, send an email to slim-discuss...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/slim-discuss/46e3e97d-eb3c-41af-8b88-108cc6abf705n%40googlegroups.com.

[Ben Haller]

Hi Sherif!

[ I see Peter has replied, and his answer has the virtue of being short.  Here's a longer answer from me, since this seems worth discussing a bit.  :-> ]

OK, so.  This is a great question, but a little bit tricky.  I'm going to take it in sort of reverse order.  :->  My answer might be more detailed than you need, but others may also find this thread useful.

First of all, if you have an x and y displacement, dx and dy, you can of course just add those in to the individual:

ind.x = ind.x + dx;
ind.y = ind.y + dy;

That works vectorized as well, if "ind" is a vector of individuals and dx and dy are vectors of per-individual displacements.  Always good to vectorize.  You can similarly use the spatialPosition property and setSpatialPosition() method to move individuals, if you prefer to work with points rather than separate coordinates.

So, where do the displacements dx and dy come from?  With a Gaussian dispersal kernel, you can draw them independently, like:

dx = rnorm(1, mean, sd);
dy = rnorm(1, mean, sd);

And again you can vectorize that by passing a count greater than 1 to rnorm().  But a Gaussian dispersal kernel is actually the *only* dispersal kernel where you get a radially symmetric result by drawing the x and y components independently; for all other dispersal kernels, you'll want to use an alternative method to get the x and y displacements.  One simple approach is to draw the angle and the radius separately, like so:

theta = runif(1, 0, 2*PI);
d = runif(1, 0, 1);

This draws an angle from 0 to 2*PI, and a distance from the uniform distribution [0,1], just for example.  It's important to note that the example above does *not* give a uniform spatial distribution within a radius of 1, however.  Here's a little R code to produce a plot illustrating that point:

thetas <- runif(10000, 0, 2*3.14159)
distances <- runif(10000)
x <- cos(thetas)*distances
y <- sin(thetas)*distances
quartz(width=3, height=3)
par(mar=c(0.5, 0.5, 0.5, 0.5), family="serif", mgp=c(2.1, 0.6, 0))
plot(x, y, pch=19, cex=0.2, xaxt="n", yaxt="n")

The distance draws are from a uniform distribution, but the resulting dispersal kernel is *not* spatially uniform, but rather clustered around the central point.  To make the distribution uniform, you can take the square root of the distances by doing:

distances <- sqrt(runif(10000))

instead.  Both of these options are perfectly valid dispersal kernels, and you might even refer to both of them as "uniform dispersal kernels", so this can be a point of confusion; be sure that you define your dispersal kernel the way you really want it to be!

OK, and now we get to what is maybe the question you really want me to answer: if you want the distance to be drawn, not from runif() or rexp() or whatever, but from your own custom distribution, how do you do that?  How can you do draws from an arbitrary distribution, in Eidos?

There are lots of answers to that, and at this level, those answers are going to be similar to the answers for the same question in R:

- Define the distribution in terms of a distribution that is supported.  For example, in R's VGAM package the rrayleigh() function, to draw from a Rayleigh distribution, is defined as:

rrayleigh <- function(n, scale = 1) {
    ans <- scale * sqrt(-2 * log(runif(n)))
    ans[scale <= 0] <- NaN
    ans
}

You could do much the same thing in Eidos, with a user-defined function or just directly in your code.

- Generate a whole bunch of draws somewhere else (R, Python, C, Mathematica...), put them in a file, read that file in at the start of your Eidos script and define it as a big vector constant, and generate draws from that dataset using sample().  This will of course only give you draws from the values in your dataset, but for many practical purposes, a dataset of 10,000 values or something will be quite an acceptable approximation, and this approach is very easy to implement.

- If you have a definition of the PDF for the distance distribution, either as a function you can code or as a binned histogram-style approximation, you can generate draws from that distribution using rejection sampling (https://en.wikipedia.org/wiki/Rejection_sampling) quite easily.

- As Peter suggests, use the inverse CDF, if you can script that.

- If the GSL (GNU Scientific Library) has a built-in function for your distribution (https://www.gnu.org/software/gsl/doc/html/randist.html), or has other math functions that Eidos doesn't have yet but that would allow you to implement your distribution, you can open a GitHub issue on SLiM to request that I add that function to Eidos.

So, there are some thoughts.  As Peter says, it's a big topic, but one that does come up often, so it's good to have a discussion of it on the list; thanks for the question.  :->

Cheers,
-B.

Benjamin C. Haller
Messer Lab
Cornell University


Sherif Negm wrote on 4/16/24 11:02 AM:

 I have a custom density function for dispersal in a continuous habitat. How can i write it down, so individuals can disperse using it. According to the manual, it only showed how individuals can disperse using preset functions like normal distribution density.