Math 3MB3 Final Project - Topic Descriptions

Each project includes a description of a base model and some possible extensions. Your group will need to first fill in the details of the base model and analyze it, then choose ONE of the suggested extensions and modify your model from base model to explore that extension.  Unless otherwise specified, you may construct either a discrete or continuous model.

Project 1:  Carbon Cycle

Carbon is a key element in terrestrial ecosystems.  It enters the soil when plants die, or shed leaves and branches (called “litterfall”, https://en.wikipedia.org/wiki/Plant   lit ter); it leaves the system by being turned into carbon dioxide by bacterial metabolism or other chemical processes. In recent history, human industrial activity has begun to interfere substantially with the biological carbon cycle in various ways, leading to potential climate change (https://en.wikipedia.org/wiki/Climate   change).

In this project, you will model carbon cycling (https://en.wikipedia.org/wiki/Carbon cycle) from the atmosphere and back through various components of an ecosystem (e.g., plants, litterfall, humus [https://en.wikipedia.org/wiki/Humus]) and examine the efects of human activity on this cycle.

Base model

Start simply by tracking carbon levels in the litter on a forest floor alone (litter being naturally-occurring debris like leaves, branches, and deadfalls, that decay over time, not human-made trash).  Let ct   be the density of carbon in the litter at time t  in months (measured in grams of carbon per square metre, or gC/m2 ).

Assume carbon enters the litter through litterfall continuously at a constant monthly rate, z, and leaves at a rate proportional to the amount of litter currently in the system via humification (the conversion of litter into humus).  Set up a simple discrete-time model for this scenario and analyze it.

Now account for atmospheric carbon by assuming that humus and litter both respirate carbon dioxide in to the atmosphere at rates proportion to the amount of each substance in the system, and that atmospheric carbon gets converted back into litter via plant growth, which produce litter at some rate proportional to the number of trees.  Set up a model for this scenario. You may find it easiest to set up one state variable for the density of carbon in each component of this system: litter, humus, atmosphere, and plants.

Finally, add some human activity to the model, and explore the repercussions.  For instance, what if humans were to plant trees each year to draw out more atmospheric carbon? What if humans were to cut trees down to use this natural resource to produce goods? How might this afect model predictions?

Possible extensions

Expanded ecosystem

Model more components of the carbon cycle: plants (subdivided into leaves, branches, stems, and roots), litter, humus, and stable humus charcoal, denoted by x1 , x2 , ...x7 , respectively. Atmospheric carbon 丑ows into plants via photosynthesis.  Leaves, branches, and stems in- crease the carbon in litter, which then increases the carbon in the humus via humification. Roots increase the carbon directly in the humus (not via litterfall).  Humus increases the carbon in the stable humus charcoal via carbonization.  Litter, humus, and stable humus charcoal all increase carbon in the atmosphere via respiration.

As a simplification, assume that the atmosphere has a constant carbon content (un- changed either by giving carbon to plants or by absorbing carbon from the litter, humus, or stable humus charcoal), since it contains so much more carbon compared to the other components of the system. You can therefore imagine that the atmosphere is outside of the system and only model the other components.  The atmosphere simply introduces a con- stant rate of carbon, z, into the system, and proportions p1 , ..., p4  indicate how much of this constant influx of carbon gets allocated to leaves, branches, stems, and roots, respectively. Parameters kij   give the rate of carbon flow from xi  to xj .  Any carbon going back into the atmosphere can be thought of as simply leaving the system, which occurs at rate ki0  for carbon flowing from xi  to the atmosphere.

A compartmental diagram of this model is as follows:

Model parameters for various ecosystems can be taken from the following table, which comes from ”A Simulation Study for the Global Carbon Cycle, Including Man’s Impact on the Biosphere” by Goudriaan and Ketner:

(The unit Gt is Gigatonnes, or a billion tonnes, where 1 tonne (metric ton) = 1000 kg =

1 Mg.)

Note that the “leaving litter” flow does not distinguish between carbon leaving litter by humification or respiration.  The humification factor h gives the proportion of the “leaving litter” flow that goes into the humus, which leaves 1-hof the flow togo into the atmosphere via respiration. Similarly, the “leaving humus” flow must be divided into the proportion that goes into the stable humus charcoal via carbonization, c, and the proportion that goes into the atmosphere, 1 - c.

You could consider one or several ecosystems based on the parameters in the above table, and compare the predicted carbon-cycling behaviour between ecosystems.

Seasonal parameters

Many parts of the world experience seasonal climate variation, which afects the growth and decay of plants.  Consider some seasonal variation in plant-related flows (e.g.  by assuming that some or all of per-capita plant-related flow rates are not constant but vary sinusoidally with a period of a year). How does this change your model predictions?

Project 2:  Drugs in the Body

The way in which drugs are administered to individuals and then metabolized by the body is of great concern in pharmacology and medicine:  drug dosage over time must be high enough for some period to have a medicinal efect on the patient, but not so high that they may overdose. There are various drug delivery methods available to those designing patient therapies. Metabolic pathways also vary depending on drug administration and type.

In this project, you will compare the efects of various drug delivery methods and metabolic pathways on a patient’s dosage over time.

Base model

The rate at which the body processes drugs depends on two factors: the rate at which the drug is administered and the drug processing rate.  Let A(t) be the amount of a drug in the body (in milligrams) at time t.  Then we can denote the dosing rate D(t), which is independent of the current amount of drug in the body A(t), and the processing rate P (A), which does depend on A. Then the model for A(t) can be cast as an “inflow minus outflow” relationship:

As a first pass, assume that a total of Dtot   mg of the drug is given intravenously at a constant rate of r mg/hour for the first h hours, at which point drug administration ceases:

We also assume that drug processing occurs linearly, proportionate to the amount of drug currently in the system: P (A) = cA, where c is the clearance rate.

Note that, while D(t) is piecewise, you can split the domain into t ∈ [0, h] and t ∈ (h, ∞), and then study what is happening in each piece of the domain to get a sense of the overall dynamics.

Assume the drug is given as a pill that is designed to dissolve slowly. As the pill dissolves, less and less of the medication is released. At any instant, the rate of release of the drug can be modelled by D(t) = Dmaxe-t/h. Calculate how much of the drug is released by the pill in the time interval [0, h].  Compare this dosage plan to the intravenous plans you previously explored.

Possible extensions

Logistic drug metabolism

Assume that instead of P (A) = cA  (exponential drug clearance), the drug is metabolized logistically, i.e., P (A) = cA(1 - A/K).  Explain how one might interpret the new logistic parameter K in this context. How does this change the amount of drug in the bloodstream over time, compared to the intravenous treatment plans explored previously?

Intermediate absorption compartment

Consider a compartmental model where the drug needs to enter another compartment before becoming bio-available (e.g., it is injected in the blood but needs todifuse into organs before it can be used, or it is swallowed into the digestive system and needs to difuse into the blood stream).  Model the amount of drug in each part of the body over time (blood, digestive system, organs—whichever apply in your context) as well as the difusion processes between these body parts.  Make a variety of assumptions for the difusion mechanisms and explore how these afect model results.

Project 3:  Conservation and Wildlife Management

It is important for ecosystems that the organisms within it remain in balance. If there are too few members of species X, it could go extinct, which may cause problems for other species that rely on X  as a food source.  In this case, humans may undertake conservation eforts to save species X . On the other hand, if species X is too abundant, it may deplete its main food source, which could pose a risk to other species that feed on the same organism. Here, humans may consider wildlife management strategies, like hunting, to keep species X from growing too abundant.

In this project, you will model an organism’s population and the efects of various con- servation and wildlife management strategies on the that population.

Base model

Many species of wildcats are endangered, including the bobcat. To better understand bobcat population dynamics, we will construct an age-structured model for a bobcat population with 16 age classes. This model should include survival and reproduction, as in the rabbit model explored in lecture. Use the following parameters in your model:

Note that, just as in lecture, we are only tracking bobcats that can give birth as a simplification, and these parameters re丑ect that assumption.

Draw a compartmental model diagram and use it to derive a general matrix model for this context (“general” meaning with symbols for parameter values) of the form(t+1) = M(t).

Use the best and worst case parameters in the model to determine the long-term model behaviour under each scenario.  What is the long-term growth rate and distribution of the population in each case?  Note there are diferent combinations of best and worst case parameters.

Now consider an intervention strategy where constant numbers of young bobcats (aged 0-2) are added to the population each year through breeding and conservation programs. How does this change your model predictions?

Possible extensions

Hunting

Assume now that instead bobcats are overly abundant and the local conservation authority is worried that they will overly consume their prey and potentially endanger their prey species, and in longer term, bobcats themselves due to extinction of their prey. Model the scenario where fully matured adult bobcats (bobcats aged 3 and older) are now instead hunted bu humans at a constant rate each year.  What does the model predict will happen in the long term?  Suppose instead that the hunting rate is proportional to the current bobcat population.  How much hunting can the bobcat population sustain before it too becomes endangered (i.e., its long-term growth rate predicts a decline)?

Catastrophes

In natural habitat of wild animals, catastrophes like wildfire or storms might happens in a regular basis and have a great impact to the population of species. Model the occurrence of sudden catastrophes to the bobcat population every n years. Assume catastrophes will not directly afect population size or survival rate of bobcats, instead they reduce the resources for bobcats, lower each age’s reproduction rate by some constant proportion p(which means the new reproduction rate is (1 - p) proportion of the original rate) just for that year. For simplicity, assume p is the same for all ages. Explore various parameter values for n and p and discuss the resilience of the population.

Project 4:  Ebola Transmission

The way in which an infectious disease spreads in a human population depends on many factors, some of which are specific to the disease and population being studied. For instance, in the 2014 West African Ebola epidemic, some disease transmission occurred at funerals, as a result of contact between susceptible individuals and individuals who had recently died of Ebola.

In this project, you will extend a simple disease model to account for transmission at funerals. You may also consider other extensions appropriate for modelling the 2014 West African Ebola epidemic.

Base model

We explored the SIR (susceptible/infected/recovered) model in lecture (we didn’t explicitly discuss the R compartment though).  There are hundreds of variants of the SIR model dealing with various complexities of disease biology and human society.  One variation is the SHERIF model http://www.sciencedirect.com/science/article/pii/S17554 36517300233, developed to analyze the recent West African Ebola outbreak, which adds Hospitalized, Exposed, and Funeral compartments to the SIR model (the order is chosen for pronounceability). Note that this paper discusses a continuous version of the model, instead of the discrete models we discussed in class.

To make things simpler, consider the SIFR model, which includes transmission caused by contact occurring at funerals.  A discrete version of this model can be encoded in the following system of equations.

St+1  = St - βSt It - βFStFt                                                                        (3a)

It+1  = It + βSt It + βFStFt - γIt                                                              (3b)

Ft+1  = Ft + γIt - σFt                                                                              (3c)

Rt+1  = Rt + σFt                                                                                      (3d)

Explain what each term in the model equations denotes, and what the parameters rep- resent.  Explore the long-term behaviour of this model.  What does it predict for various parameter values?

Add vital dynamics, i.e., births and deaths, to the model: St+1  gains a +μN - μSt  term; the other compartments use a loss term with a per capita rate μ (e.g.  -μIt   for the It+1 compartment)) and explore the long-term behaviour of the model.  (If necessary, you can drop the F compartment and analyze the SIR equations.)

Possible extensions

Controlling funeral spread

Model public health interventions that reduce transmission at funerals by diverting some people from I directly to R without going through F. How does this change model predic- tions?

Add exposed and hospitalized individuals

Add the E and H terms to produce a discrete version of the full SHERIF model. Compare this model’s results to that of the SIFR model analyzed previously (with and without vital dynamics).

Project 5:  Population predictions

In order to plan for the future, governments often use projections of population size to project demands on social services.  Various factors afect population size: births, deaths, immigration and emigration. In this project, you will build a demographic model to project Ontario’s population growth or decay in the future, based on population models of California, US.

Base model

The base model for this project will focus on population changes in California compared to the rest of the United States. You will then extend this to do something similar for Ontario’s population compared to the rest of Canada, and compare between Ontario and California.

Let (t) denote the populations of California and the US excluding California at time t, respectively (so the first component is California’s population and the second is the rest of the US). The following table gives population values in 1955 and 1960, along with births, deaths, and net migrations (immigration - emigration) between these two regions: 

Note that these population counts are in units of 1,000 people.  Use this information to calculate the per-capita rates associated with each process by filling in the following table: 

Two rates have already been computed for you so that you can verify your calculations are sound. Set up a matrix model of the form. (t + 1) = (I + B - D + M)(t),where the matrices B, D, M are birth, death, and net migration matrices, respectively (which should all be diagonal), and I is the identity matrix. Explain why I is needed here.

Using the 1955 data and the rates you computed, verify that your model accurately pre- dicts the 1960 population.  Simulate the model to 2020 and verify your predictions against the https://www.census.gov/data.html.  How close were you predictions to the ac- tual outcomes? What may have caused discrepancies between you projection and the true outcomes?  What does your model predict in the long run for the population distribution between California and the rest of the US (as t → +∞)? You may find it easier to perform. this analysis if you define a general ”growth” matrix G  = I + B - D + M.  What is the long-term rate of growth of the population of California? What proportion of US residents will live in California in the long run, according to your model?

Now work with immigration and emigration separately (instead of looking at net migra- tion). Migration between the two regions is summarized in the following table:

Ignoring births and deaths, construct a matrix model of transitions between regions

with the form. (t + 1) = T(t).  Again, simulate the model to 2020 and compare it with

US census estimates.  Discuss any discrepancies.  Determine the long-term behaviour of the model, including the long-run growth rate of California’s population, as well as the population distribution within the US. How do these results compare to the earlier model with births, deaths, and net migration?

Possible extensions

Comparing to Ontario and Canada

Apply your model to Ontario and the rest of Canada.  Look for realistic parameter values using census data from Statistics Canada https://www12.statcan.gc.ca/census-recen sement/index-eng.cfm. Use your model to make projections for the population past 2022. How do your estimates compare with those from Statistics Canada https://www150.statc an.gc.ca/n1/en/catalogue/91-520-X? Explore their projection methodology to explain any diferences between your projections and theirs.