bodkan / slendr

Population genetic simulations in R 🌍
https://bodkan.net/slendr
Other
54 stars 5 forks source link

fixes for SLiM 4.1 changes to SpatialMap #147

Closed bhaller closed 9 months ago

bhaller commented 9 months ago

Hi @bodkan! I figured I might as well open an issue for this:

WARNING: SLiM 4.1 which has just been released includes a couple of backwards incompatible changes related to the implementation of spatial maps which prevent the current version of slendr's slim() function from working correctly. If you rely on the functionality provided by the slim() function, you will have to use SLiM 4.0. (Note that if you want to have multiple versions of SLiM on your system, you can either use the slim_path = argument of slim() or specify the $PATH to the required version of SLiM in your ~/.Renviron file just like you do under normal circumstances). Porting slendr for SLiM 4.1 is being worked on.

Mostly I'm just curious: what did I break? As far as I recall now, my changes to SpatialMap were all backwards-compatible; I certainly intended them to be. As the release notes for 4.1 say:

It should preserve backward compatibility with SLiM 4.0 and 4.0.1 (models that run in those will continue to run).

Sorry if I broke something inadvertently!

bodkan commented 9 months ago

Ughh, so this is a little embarrassing -- right now I don't remember what exactly broke!

What happened is that when I upgraded to SLiM 4.1 (right after you released it) and re-ran my unit tests, the SLiM back-end script inside slendr started crashing. I just remember the crashes/errors pointed to something related to the spatial maps... but not more.

I decided to deal with this at a later point after this release (it's nothing dramatic) but then I had to postpone the release because viruses & travel... and I got to it just now, apparently having forgotten in the meantime what exactly did break.

bodkan commented 9 months ago

Sorry if I broke something inadvertently!

And no need to apologize! There's a reason why I didn't ping you before I do my own investigation. :) I haven't had a chance yet to carefully read your changelog. I just thought "OK, some minor updates to the spatial maps; will investigate later".

But you know what, since you've opened the issue (thank you!), let me run something with slim() under SLiM 4.1 right now and report back. I'd have to do this anyway, so at least I'll have something to think about. And I'm curious myself.

bodkan commented 9 months ago

OK, here's a minimal little R script that crashes under SLiM 4.1 but runs under SLiM 4.0.1. You don't have to read that, @bhaller, it's mostly for my own reference.

library(slendr)
init_env()

# create an example map
map <- world(xrange = c(-25, 55), yrange = c(-32, 35), crs = 4326)

# create a single population spread all over the map
pop <- population("pop", time = 1, N = 2000, map = map)

# compile a single-population model
model <- compile_model(
  populations = pop, generation_time = 30, simulation_length = 1000,
  resolution = 1, competition = 1, dispersal = 1, mating = 1
)

# run the model through SLiM, return a tree sequence
ts <- slim(model, sequence_length = 1000, recombination_rate = 0)

Here's a full long from the SLiM run which answers your question:

// Initial random seed:                                                                          
7837909936800388019

// RunInitializeCallbacks():
SEED: 7837909936800388019
initializeSLiMOptions(keepPedigrees = T, dimensionality = 'xy');
initializeInteractionType(0, "xy", reciprocal=T, maxDistance=1);
initializeInteractionType(1, "xy", reciprocal=T, maxDistance=1);
initializeTreeSeq();
initializeMutationType(0, 0.5, "f", 0);
initializeGenomicElementType(1, m0, 1);
initializeGenomicElement(g1, 0, 999);
initializeMutationRate(0);
initializeRecombinationRate(0);

// Starting run at tick <start>:
1 

Generation 1: starting the simulation
Generation 1: creating pop(p0)
Generation 1: updating map for pop(p0)
ERROR (Subpopulation::ExecuteMethod_setSpatialBounds): setSpatialBounds() new spatial bounds are not compatible with an attached map named 'world'; use removeSpatialMap() to remove incompatible spatial maps before changing the spatial bounds.  (This enforces internal consistency and avoids accidentally stretching a map to new spatial bounds.)

Error on script line 222, character 12:

        pop.setSpatialBounds(c(0.0, 0.0, asFloat(WIDTH) - 1, asFloat(HEIGHT) - 1));
            ^^^^^^^^^^^^^^^^
Error: Unfortunately SLiM terminated before a tree sequence was saved.
See the above for an indication of where things could have gone wrong.

The error happens on this line of the SLiM script which is bundled with slendr. As far as I can remember (not seeing that code for literal months 😬), the surrounding s4 callback is called to either define (upon creating a population for the first time) or update (whenever a spatial population boundary changes) a binary raster spatial map of a population.

Interestingly enough, I don't remember seeing this particular error. I thought it was the pop.defineSpatialMap() call just above that line 222 that was giving problems... then again, as I wrote above, it's been weeks since I did that trial SLiM 4.1 test run and must've forgotten what was the problem.

Also, the unit tests which are now crashing are actually all crashing on that same line 222 (so I did clearly misremember what was going on). So, it's clearly not an issue with any backwards incompatible change. Just SLiM being more strict to avoid user errors.

Now, I don't exactly remember why am I calling setSpatialBounds() after each map definition with defineSpatialMap(). Mind you, these lines of code represent the most ancient bits of slendr before it was even slendr. In fact, before there was even any idea of putting together a spatial simulation R package. :) This was possibly written just hours after I first read in the SLiM manual that SLiM can do simulations in a continuous space with a raster map. :D

Am I wrong or should the setSpatialBounds() be called only once as soon as a population is created, but it's not necessary (in fact, it's now wrong) to call it again after spatial raster maps are updated? Looking at this code, it seems to me that I thought that when a population map is defined, the user should also set the boundaries... but that's probably not true.

bodkan commented 9 months ago

Note (mostly) to self: If my intuition is correct, it will be trivial to fix this by moving the setSpatialBounds() call somewhere here and that's it. And, of course, correcting the changelog, to avoid spreading unsubstantiated panic.

bhaller commented 9 months ago

Hi @bodkan. Aha, OK. This is a break in backward compatibility, yes, but it is regarded as more or less a bug fix. The error indicates something that probably should always have been an error, but previously wasn't, and now is. This error occurs when you have defined a spatial map on a subpopulation, and then you change the bounds of the subpopulation. Spatial maps consider themselves to be associated with a specific set of bounds, and require that all subpopulations that they are set on match those bounds exactly. This keeps the spatial map and all subpops using it tied to the same coordinate system, which is essential to the coherence/integrity of the simulation. If these requirements are too stringent in some way, and need to be relaxed somewhat for your purposes, I could consider that. (Perhaps the origin of the bounds could differ, as long as the size matches?) But the overall policy seems correct to me. Open to discussion.

bhaller commented 9 months ago

(If you set the bounds again to the same bounds you used before, though, you should not get an error. This error indicates that you are changing the bounds after the spatial map was created and associated with the subpop; that is what is not allowed.)

bodkan commented 9 months ago

Hi @bodkan. Aha, OK. This is a break in backward compatibility, yes, but it is regarded as more or less a bug fix. The error indicates something that probably should always have been an error, but previously wasn't, and now is. [...] If these requirements are too stringent in some way, and need to be relaxed somewhat for your purposes, I could consider that. (Perhaps the origin of the bounds could differ, as long as the size matches?) But the overall policy seems correct to me. Open to discussion.

OK, that's what I assumed is happening and I agree with your motivation. And no, I don't think this is a restriction that needs to be relaxed to a previous state just because of this -- I'm not relying on this behavior at all.

Let me just try moving the setSpatialBounds() call to just after each population is created (not after each map update). If I understand everything correctly, this should fix the error without changing anything in terms of model behavior.

With that being said:

(If you set the bounds again to the same bounds you used before, though, you should not get an error. This error indicates that you are changing the bounds after the spatial map was created and associated with the subpop; that is what is not allowed.)

In the call pop.setSpatialBounds(c(0.0, 0.0, asFloat(WIDTH) - 1, asFloat(HEIGHT) - 1)); which happens here, note that WIDTH and HEIGHT are constants set in the initialization() block here. In other words, these are fixed based on the dimension of the very first spatial raster passed from slendr and the very beginning of the SLiM run (and all such rasters are of fixed dimension).

So, the bounds being re-set repeatedly (which causes the error) are, in fact, the same bounds used again and again. So if you don't consider it a problem if setSpatialBounds() is called with the same bounds then perhaps SLiM 4.1 is a little stricter than you intended?

Then again, as I said, this is not a problem for slendr and not something I rely on.

bhaller commented 9 months ago

OK, that's what I assumed is happening and I agree with your motivation. And no, I don't think this is a restriction that needs to be relaxed to a previous state just because of this -- I'm not relying on this behavior at all.

Let me just try moving the setSpatialBounds() call to just after each population is created (not after each map update). If I understand everything correctly, this should fix the error without changing anything in terms of model behavior.

OK, sounds good.

With that being said:

(If you set the bounds again to the same bounds you used before, though, you should not get an error. This error indicates that you are changing the bounds after the spatial map was created and associated with the subpop; that is what is not allowed.)

In the call pop.setSpatialBounds(c(0.0, 0.0, asFloat(WIDTH) - 1, asFloat(HEIGHT) - 1)); which happens here, note that WIDTH and HEIGHT are constants set in the initialization() block here. In other words, these are fixed based on the dimension of the very first spatial raster passed from slendr and the very beginning of the SLiM run (and all such rasters are of fixed dimension).

So, the bounds being re-set repeatedly (which causes the error) are, in fact, the same bounds used again and again. So if you don't consider it a problem if setSpatialBounds() is called with the same bounds then perhaps SLiM 4.1 is a little stricter than you intended?

I don't think there's a problem here; I just checked. Setting the bounds repeatedly to the same thing does not trigger the error message. Either you're setting the bounds inconsistently, or – more likely, I'm guessing – when you first create the spatial map the bounds are still the default (0, 0, 1, 1) bounds from SLiM, and so when you later set the bounds to something different, they are inconsistent with those original (0, 0, 1, 1) bounds. The bounds need to be correct already when you create the spatial map, because the spatial map inherits the bounds from the subpopulation at that point.

bodkan commented 9 months ago

[...] when you first create the spatial map the bounds are still the default (0, 0, 1, 1) bounds from SLiM, and so when you later set the bounds to something different, they are inconsistent with those original (0, 0, 1, 1) bounds. The bounds need to be correct already when you create the spatial map, because the spatial map inherits the bounds from the subpopulation at that point.

Ah, yes. That makes sense and must be, given how the SLiM script in slendr is written, exactly what's going on.

The unit tests have passed and I'm re-rendering all examples on the website just to make sure there's no other issue.

Thanks for your help!