# Resource-worm-immune model

[https://doi.org/10.5061/dryad.59zw3r2gx](https://doi.org/10.5061/dryad.59zw3r2gx)

## Overview

Users will usually:

1.  Set parameters in RWI_Parameters.R
2.  Use RWI_Dynamics.R to calculate trajectories for various parameters sets and plot them as heat maps (time on x-axis, host investment on y-axis, dependent parameter on z-axis)
3.  Use RWI_Compare.R to calculate optimal host investment in different strategies under different parameter values
4.  Use RWI_PlotOutput.R to plot results of using RWI_Compare.R

The code uses R packages tidyverse (Wickham et al. 2019) and deSolve (Soetaert et al. 2010).

## Code files (uploaded to Zenodo)

##### RWI\_Parameter.R

Contains default parameter values and mathematical models needed to perform simulations. See the associated article for a model description and derivation.

*   The time period over which to carry out simulations can be changed by resetting tmax.
*   For some parameters, solutions are required over a range of parameters. Specifically: host investment c, resource availability SR, infection pressure SL, larval and adult parasite harm hL and hM (currently set so hL+hM=0.8), parasite maturation rate g and immunopathology hI. The name for each range is given by appending "s" to the associated parameter name.
*   Default parameter values are also given. Commented out values are those used for parameter sensitivity tests.
*   f_check_w checks that parameters are sensible by calculating time to host death if starved; this value is returned in the command line. Also returned in the command line are quantities which must be positive for the host to initially be alive, and the associated values of SR. These checks allow the user to ensure that the parameters chosen are at least plausible. See Appendix of associated article for further details.
*   f.sys is the ODE system for single strategies. Users are unlikely to need to directly use this. 
*   f.traj uses f.sys to calculate the model trajectory. Users must specify strat (host parasite-mitigation strategy), parms.f (named vector of parameters differ from default values; can be NULL) and t.spec (vector of time points at which to return results).
*   f.sys.combine is the ODE system for a combined strategy. Users are unlikely to need to directly use this. 
*   f.traj.combine uses f.sys.combine to calculate the model trajectory. Users must specify parms.f (named vector of parameters that differ from default values; can be NULL) and t.spec (vector of time points at which to return results).
*   Note that a single strategy is represented by setting the variable "strat" to 1 (anorexia), 2 (target larvae), 3 (target established parasites), 4 (tolerance).  "strat" can also be set to 0 (no strategy) or 1.5 (starvation). When modelling a combined strategy, instead the allocation of immune response between strategies must be specified (or optimised for).

##### RWI\_Dynamics.R&#x20;

Code allowing the user to calculate model trajectories for each single strategy over a range of values of host investment c, infection pressure SL and resource availability SR. This can easily be adapted to consider different parameter ranges, as the trajectory functions in RWI_Parameters.R allow the user to pass whichever parameters they desire to specify via a named vector.

*   sol.x (for x=i,ii,iii,iv) are dataframes containing the trajectories for each strategy (1,2,3,4).
*   f.DeadHost replaces system values with NA when the host has dead (condition reached zero).
*   The code loops through each parameter value and calculates the trajectory for each strategy.
*   The dataframes of solutions are then normalised and put into "long" form to facilitate plotting.
*   c.opt is the output from a previously computed optimisation. This is included in the code in order to indicate optimum host investment on heat maps in which the x-axis is time, y-axis is host investment and z-axis is the value of the appropriate dependent variable.
*   sol.plot.t.all.X is a function which plots the results for a particular SR and SL, for a choice of dependent variable (input X; taking values C,R,L,M, or I, representing condition, host resources, larval burden, established parasite burden, or immune response).

##### RWI\_Compare.R&#x20;

Contains code allowing the user to optimise host investment in immune strategies over a range of parameter values. This uses the code contained in RWI_Optimise.R and RWI_Combine.R.

*   Running RWI_Optimise.R outputs c.opt, a dataframe containing results. 
*   Sometimes particular parameter combinations fail to converge (c.opt$convergence != 0); a process for finding these and retrying with a (manually entered) different initial guess is coded here.
*   Sometimes the optimisation converges at a local rather than global maximum. These points must be found by eye; a process for finding these and retrying with a (manually entered) different initial guess is coded at the end of this file.
*   f.MinSurvival calculates the minimum value of SR at which hosts survive, and adds data allowing this to be plotted (see Figure 5 in the associated article, for example).
*   The same process is also carried out for optimising the combined strategy, using RWI_Combine.R. In this case the host investment and the allocation between the two resistance and tolerance strategies are optimised concurrently.

##### RWI\_Optimise.R&#x20;

Calculates the optimum host investment in each strategy, over the specified resource availabilities SR, infection pressures SL, parasite maturation rates g and larval harms hL. Note established parasite harm hM is constrained so that hL+hM=0.8.

*   trap.rule is a function that estimates an integral over a vector using the trapezium rule.
*   f.optim.fn is the function to optimise over.
*   f.optim carries out the optimisation. Users must specify the strat (host strategy), parms.f (named vector of parameters that differ from default values; can be NULL) and initial guess for the optimisation (c0). The function returns a dataframe containing strategy, value of c, mean host condition C.Mean, final values of R,M,L,I,C, and an indicator of convergence (0 if converged).

##### RWI\_Combine.R

Same as RWI_Optimise.R, but for the combined strategy.

##### RWI\_PlotOutput.R

Produces plots such as Figure 5 in the associated article.

*   Uses f_plot_opt.R to plot dependent parameters from the output of RWI_Optimise.R. Users must specify data to plot, (dat; output from RWI_Optimise.R), dependent parameter to plot (Y), infection pressure (SL; one of the values in the SL column in dat), immunopathological harm (hIs; one of the values in the hI column in dat), title (default NULL) , XLim (x-axis limits; default c(0,6)), YLim (y-axis limits, default c(0,4)), y-label (YLab; default NULL), combined dataset to include (output from RWI_combine.R; must have same parameter values as dat, default NULL)
*   Uses f_plot_opt_nu.R to plot the allocation of the immune response in the combined strategy

## Simulation data files (uploaded to Dryad)

Data files contain simulation results used in the associated article. These are outputs of RWI_Optimise.R and RWI_Combine.R, calculating optimum host investment in each strategy over a range of parameter values. 

For hosts, with a single strategy, columns are labelled:

*   strat: the strategy used, taking values 1 (anorexia), 2 (larval prevention), 3 (adult clearance), 4 (tolerance)
*   c: the optimum host investment for this set of parameter values
*   SR: resource availability
*   SL: infection pressure
*   g: parasite maturation
*   hL: larval harm
*   hM: adult parasite harm, with hL+hM=0.8
*   hI0: linear coefficient of immunopathological harm
*   hI1: quadratic coefficient of immunopathological harm
*   C.Mean: mean host condition over simulation time period
*   R: host resources at end of simulation
*   M: adult parasite load at end of time period
*   L: larval parasite load at end of time period
*   I: immune levels at end of time period
*   C: host condition at end of time period
*   convergence: indicator of convergence (0 if converged; all solutions in these datasets have converged)

For hosts with a combined strategy, there is no "strat" column. Instead there are the following:

*   nuL: proportion of immune response allocated to preventing larval establishment
*   nuP: proportion allocated to clearing established parasites
*   nuT: proportion allocated to tolerance (ameliorating parasite harm) 

Unless otherwise indicated, the simulations took place over 90 days.

Files are as follows:

*   opt_harm0p8.csv contains results for the default parameter values (found in RWI_Parameters.R)
*   opt_harm0p8_q0p5.csv, opt_harm0p8_v2p5.csv, opt_harm0p8_v2p5_q0p5.csv. contain the same, but with (respectively) immune production q=0.5, immune saturation coefficient v=2.5, and both q=0.5 and v=2.5.
*   opt_harm0p8_k0_5.csv, opt_harm0p8_k1_2p5.csv, opt_harm0p8_k0_5_k1_2p5.csv. contain the same, but with (respectively) immune strength  k0=0.5, immune strength saturation coefficient k1=2.5, and both k0=0.5 and k1=2.5.
*   opt_harm0p8_a4w2.csv, opt_harm0p8_r5.csv, opt_harm0p8_a4w2_r5.csv. contain the same, but with (respectively) condition conversion a=4 and condition loss w=2, resource consumption r=5, and all of a=4, w=2 and r=5.

<!---->

*   opt_harm0p8_Combine.csv contains the optimisation output for the combined strategy under default parameter values.
*   opt_harm0p8_t7.csv contains the optimisation output for single strategies under default parameter values over 7 days.

## References

Karline Soetaert, Thomas Petzoldt, R. Woodrow Setzer (2010). Solving Differential Equations in R: Package deSolve. Journal of Statistical Software, 33(9), 1--25. doi:10.18637/jss.v033.i09

Wickham H, Averick M, Bryan J, Chang W, McGowan LD, François R, Grolemund G, Hayes A, Henry L, Hester J, Kuhn M, Pedersen TL, Miller E, Bache SM, Müller K, Ooms J, Robinson D, Seidel DP, Spinu V, Takahashi K, Vaughan D, Wilke C, Woo K, Yutani H (2019). “Welcome to the tidyverse.” *Journal of Open Source Software*, *4*(43), 1686. doi:10.21105/joss.01686 [https://doi.org/10.21105/joss.01686](https://doi.org/10.21105/joss.01686).
