Merging stochasticity functioning into develop

[Thomalpas]

Hey all!

Ok so these are I think all the changes necessary to add my stochasticity functioning into develop. I'll detail them below but I think apart from small changes to ModelParameters & simulate everything is a new function so hopefully merging shouldn't be too much of a headache... Though the changed files below show for instance using DiffEqCallbacks as an addition of mine into BEFWM2.jl but that's already there in develop (but not main), so maybe I'm doing something wrong.

I'm not great at github so if there is a way of creating a copy of the develop branch and then adding my changes to that then starting a pull request then let me know! As it is this a copy of main that I have first tried to make identical to develop and then added my changes on top of.

As a quick summary I've built a module that selects species to be stochastic based on user input, and then gives a parameter for these species (growth for producers, metabolism for consumers) its own differential equation, which along with the biomass equations are passed into an SDEProblem and solved. I'll detail more below if it helps with merging

Key changes that influence code used by everyone:

I've created a lot of new functions that are only required for stochasticty, but the below files are the only places where I think I have modified code already found in the develop branch.

• BEFWM2.jl - 3 new dependencies, 2 new files, export of 2 new functions
• model/modelparameters.jl - Incorporation of AddStochasticity module into ModelParameters function and display
• model/simulate.jl - An if loop within the simulate function that uses the unedited dbdt! - ODEProblem route as a default (if a boolean within the new AddStochasticity module is false), but switches to my new SDE solving pathway if true

Key changes that should only affect stochasticity

2 new files

• inputs/stochasticity.jl to produce the AddStochasticity module for ModelParameters
• model/covariance_matrix.jl to produce a covariance matrix that correctly supplies the CorrelatedWienerProcess providing noise with environmental and demographic stochasticity

Stochasticity-compatible alternate functions

I've modified a lot of files but mostly this has been creating new versions of functions. These new functions are only called within simulate if the boolean in AddStochasticty is true and so they shouldn't affect anything else. All the original functions mentioned below have been left untouched.
In most cases these new versions are necessary because we're now solving additional differential equations. The vector of biomass B has been expanded to hold stochastic parameters r or x and renamed diffparams (though not too important this probably needs a better name). Here are the main examples to give an idea:

• stochequations! instead of dBdt! - produces BEFW equations for the biomasses and OU processes for the stochastic parameters
• _stoch_logistic_growth instead of logisticgrowth - The r term for stochastic producers is held in diffparams not BioRates
• _stoch_consumption instead of consumption - The x term for stochastic consumers is held in diffparams not BioRates
• _stoch_metabolic_loss instead of metabolic_loss - The x term for stochastic consumers is held in diffparams not BioRates

Docs and testing files still to come but this should help give an idea of where all my stuff fits into the rest of the package!

Very happy to hear questions and possible changes 😀

[Thomalpas]

2 questions probably best for @ilajaait & @iago-lito :

1. I don't think the PositiveDomain callback will work with the CorrelatedWienerProcess that I'm using to supply the noise as the callback has to be adaptive and the noise process can't be. What do you think would be the problem with not having this callback?
2. What's the reasoning for taking the ModelParameters out of the simulate output? I've got a few functions regarding calculating dominant eigenvalues in stability.jl (which may not be suitable for including in this package anyway as it's pretty niche). These would be 'safer' if ModelParameters was in simulate output as this would ensure that the matching equilibrium biomass and model parameters are used together, which is currently not necessarily the case

[iago-lito]

Hi @Thomalpas, and congratulations for pushing all this work forward into a PR!

Regarding the merge, I can have a look into your modifications ASAP and check whether it can be done easily. From what I understand, you have based your work on main, then attempted to mimick the evolution between main and develop, and then added your own work on top of that, right? I do expect conflicts, but since you have taken the time to carefully list and explain all your changes, it'll help me a lot when navigating your work :)

Regarding your 1. concern, maybe @ilajaait can answer better than I?

Regarding 2., I'm sorry I'm not sure I understand :( From what I understand, ModelParameters values are given as input to the simulate() function, and they are not modified within it. As a consequence, they need not be returned as an output because anything simple like:

mp = ModelParameters(...)
trajectories = simulate(mp, ...)
safe_bundle = (mp, trajectories)

would ensure that your results and your parameters are "used together". Does this make sense in the context of something you meant?

[ismael-lajaaiti]

Hey @Thomalpas congrats for the PR 🎉

Concerning your issues:

1. I still don't know at 100% how behaves PositiveDomain callback (have to do some tests), but you're probably right the stochastic process and the PositiveDomain callback won't go well together (but has to be tested). Just to be sure to understand well, could you detail a little bit how the stochastic process works/is implemented? (Or how the noise term is written in the dynamical equations? You can just send the link of a paper that use a similar model). Otherwise, you can totally deactivate the PositiveDomain callback it is just a safety net to ensures that no species has a negative biomass (and adapt the simulation timesteps to do so IIUC). An option could be to create a NegativeBiomass callback whose goal would be to throws a warning (and no more) when a species has a negative biomass, so we still have a net that warns the user that something may have gone wrong during the simulation.
2. I agree with @iago-lito 's answer, but may be you mean that ModelParameters has to be modified during simulate e.g. remove extinct species (I assume that can be useful to measure the stability), no?

A note concerning the variable naming, we have to decide whether or not we put a _ internal function names. I didn't until now but we can change, doesn't matter to me, we just have to be consistent IMO. Regarding the stochastic alternative functions, may be the naming can also be more consistent e.g. original_function always becomes stoch_original_function (so stoch_dBdt! instead of stochequations! for instance).

Otherwise I'm curious about the stability metrics you have implemented because I intend to write a companion package to BEFWM2 to measure the stability of an ecological system with various metrics. 👀

[iago-lito]

Hi again @Thomalpas. I have solved conflicts and squashed all your branch into a single commit on top of develop, so it should be up-to-date now. You can force-pull stochasticity_v2 to check it on your side.

Your original commit history based on main is not deleted. I have just moved it to branch stochasticity_main_based in case we need to fetch anything back inside.

Unfortunately, there are no tests files for stochastic simulations yet. As a consequence, although I can ensure that git flow is restored, I cannot ensure that everything is working as you intend :/ Can you check on your side that everything is alright with the new squashed commit? (372c070)

[iago-lito]

Regarding internal function naming, prefixing "private" identifiers with _ is useful in e.g. python because all identifiers are exported from every module and visible to users. So users know they are doing something expected by developers when referring to a_function(), and something unexpected by developers when referring to _a_function().

In Julia, this distinction is enforced by the export keyword in the modules system. In the following situation for example:

# src/BEFWM2.jl
function public_function() ... end
function private_function() ... end
export public_function

The users are expected to tell the difference this way:

using BEFWM2 # Only bring 'export'ed identifiers into scope
public_function() # If I can use it, then it's public.
private_function() # JULIA ERROR: I can't use it, then it's private.

# But if I still want to access it, I can.
using BEFWM2: private_function # Explicit import.
private_function() # I can use it, but I know it's private.

As a consequence, the _ prefix is unnecessary imo.

[Thomalpas]

Regarding 2., I'm sorry I'm not sure I understand :( From what I understand, ModelParameters values are given as input to the simulate() function, and they are not modified within it. As a consequence, they need not be returned as an output because anything simple like:

mp = ModelParameters(...)
trajectories = simulate(mp, ...)
safe_bundle = (mp, trajectories)

would ensure that your results and your parameters are "used together". Does this make sense in the context of something you meant?

Regarding 2, that is exactly what I meant, sorry for not being clear, I wasn't sure if 'safe_bundle' as you've called it was standard procedure but it was what I was using

[Thomalpas]

Sorry I'm still getting used to GitHub

Would I want to accept all changes from left? This is your updated commit isn't it

[Thomalpas]

@iago-lito @ilajaait want to have a quick meeting in the next few days? It would probably be easier to describe the stochastic process and Ismaël to discuss stability metrics aloud rather than typing? I'm free anytime Thursday and anytime except 12:30-1pm Montpellier time tomorrow

[Thomalpas]

Actually no you're an hour ahead it's 2:30-3pm I can't do tomorrow (1:30-2 UK time)

[iago-lito]

@Thomalpas Sure, Thursday is fine for me, details on the Slack channel ;)

[iago-lito]

Status: this PR has been automatically, and definitely, closed by github because the underlying branch has been renamed to stochasticity, but the proposition remains.

[iago-lito]

Dismissing in favour of #134.