Mathematics In Population Biology Thieme Pdf
- 9,314,282 books books
- 84,837,646 articles articles
- ZLibrary Home
- Home
Main Mathematics in Population Biology
Mathematics in Population Biology
Horst R. Thieme
How much do you like this book?
What's the quality of the file?
Download the book for quality assessment
What's the quality of the downloaded files?
The formulation, analysis, and re-evaluation of mathematical models in population biology has become a valuable source of insight to mathematicians and biologists alike. This book presents an overview and selected sample of these results and ideas, organized by biological theme rather than mathematical concept, with an emphasis on helping the reader develop appropriate modeling skills through use of well-chosen and varied examples.
Part I starts with unstructured single species population models, particularly in the framework of continuous time models, then adding the most rudimentary stage structure with variable stage duration. The theme of stage structure in an age-dependent context is developed in Part II, covering demographic concepts, such as life expectation and variance of life length, and their dynamic consequences. In Part III, the author considers the dynamic interplay of host and parasite populations, i.e., the epidemics and endemics of infectious diseases. The theme of stage structure continues here in the analysis of different stages of infection and of age-structure that is instrumental in optimizing vaccination strategies.
Each section concludes with exercises, some with solutions, and suggestions for further study. The level of mathematics is relatively modest; a "toolbox" provides a summary of required results in differential equations, integration, and integral equations. In addition, a selection of Maple worksheets is provided.
The book provides an authoritative tour through a dazzling ensemble of topics and is both an ideal introduction to the subject and reference for researchers.
Publisher:
Princeton University Press
Series:
Princeton Series in Theoretical and Computational Biology, 12
The file will be sent to your email address. It may take up to 1-5 minutes before you receive it.
The file will be sent to your Kindle account. It may takes up to 1-5 minutes before you received it.
Please note: you need to verify every book you want to send to your Kindle. Check your mailbox for the verification email from Amazon Kindle.
Most frequently terms
Mathematics in Population Biology PRINCETON SERIES IN THEORETICAL AND COMPUTATIONAL BIOLOGY Simon A. Levin, Series Editor Mathematics in Population Biology Horst R. Thieme P R I N C E T O N U N I V E R S I TY P R E S S P R I N C ET O N A N D O X F O R D c 2003 by Princeton University Press Copyright Published by Princeton University Press, 41 William Street, Princeton, New Jersey 08540 In the United Kingdom: Princeton University Press, 3 Market Place, Woodstock, Oxfordshire OX20 1SY All rights reserved Library of Congress Cataloguing-in-Publication Data Thieme, Horst R., 1948– Mathematics in population biology / Horst R. Thieme p. cm. -- (Princeton series in theoretical and computational biology) Includes bibliographical references and index. ISBN 0-691-09290-7 (cl : alk. paper) -- ISBN 0-691-09291-5 (pbk. : alk. paper) 1. Population biology--Mathematical models. 2. Communicable diseases--Mathematical models. I. Title. II. Series. QH352 .T45 2003 577.8 8 015118--dc21 British Library Cataloguing-in-Publication Data A catalogue record for this book is available from the British Library. This book has been composed in Times. Typeset by T&T Productions Ltd, London. ∞ Printed on acid-free paper. www.pupress.princeton.edu Printed in the United States of America 10 9 8 7 6 5 4 3 2 1 To Adelheid Contents xiii Preface Chapter 1. Some General Remarks on Mathematical Modeling Bibliographic Remarks 1 3 PART 1. BASIC POPULATION GROWTH MODELS 5 Chapter 2. 7 2.1 2.2 2.3 The Fundamental Balance Equation of Population Dynamics Birth Date Dependent Life Expectancies The Probability of Lifetime Emigration Chapter 3. 3.1 3.2 Birth, Death, and Migration Unconstrained Population Growth for Single Species Closed Populations 3.1.1 The Average Intrinsic Growth Rate for Periodic Environments 3.1.2 The Average Intrinsic Growth Rate for Nonperiodic Environments Open Populations 3.2.1 Nonzero Average Intrinsic Growth Rate 3.2.2 Zero Average Intrinsic Growth Rate Chapter 4. Von Bertalanffy Grow; th of Body Size Chapter 5. 5.1 5.2 5.3 5.4 5.5 The Bernoulli and the Verhulst Equations The Beverton–Holt and Smith Differential Equation 5.2.1 Derivation from a Resource–Consumer Model 5.2.2 Derivation from Cannibalism of Juveniles by Adults The Ricker Differential Equation The Gompertz Equation A First Comparison of the Various Equations Chapter 6. 6.1 6.2 Classic Models of Density-Dependent Population Growth for Single Species Sigmoid Growth General Conditions for Sigmoid Growth Fitting Sigmoid Population Data 7 9 11 13 13 14 17 19 21 28 33 37 37 39 40 42 45 47 47 51 52 57 viii CONTENTS Chapter 7. The Allee Effect 7.1 7.2 7.3 First Model Derivation: Search for a Mate Second Model Derivation: Impact of a Satiating Generalist Predator Model Analysis Chapter 8. Chapter 9. 9.1 9.2 9.3 9.4 10.5 10.6 10.7 75 Discrete-Time Single-Species Models 81 The Discrete Analog of the Verhulst (Logistic) and the Bernoulli Equation: the Beverton–Holt Difference Equation and Its Generalization The Ricker Difference Equation Some Analytic Results for Scalar Difference Equations Some Remarks Concerning the Quadratic Difference Equation Bibliographic Remarks Modeling Toxicant and Population Dynamics Open Loop Toxicant Input Feedback Loop Toxicant Input Extinction and Persistence Equilibria and a Threshold Condition for Population Extinction Stability of Equilibria and Global Behavior of Solutions Multiple Extinction Equilibria, Bistability and Periodic Oscillations Linear Dose Response Bibliographic Remarks Chapter 11. Population Growth Under Basic Stage Structure 11.1 11.2 11.3 11.4 11.5 11.6 11.7 11.8 11.9 65 67 69 Nonautonomous Population Growth: Asymptotic Equality of Population Sizes Chapter 10. Dynamics of an Aquatic Population Interacting with a Polluted Environment 10.1 10.2 10.3 10.4 65 A Most Basic Stage-Structured Model Well-Posedness and Dissipativity Equilibria and Reproduction Ratios Basic Reproduction Ratios and Threshold Conditions for Extinction versus Persistence Weakly Density-Dependent Stage-Transition Rates and Global Stability of Nontrivial Equilibria The Number and Nature of Possible Multiple Nontrivial Equilibria Strongly Density-Dependent Stage-Transition Rates and Periodic Oscillations An Example for Multiple Periodic Orbits and Both Supercritical and Subcritical Hopf Bifurcation Multiple Interior Equilibria, Bistability, and Many Bifurcations for Pure Intrastage Competition Bibliographic Remarks 81 83 84 99 104 107 108 114 117 120 125 135 139 149 151 151 153 155 156 157 160 162 166 168 181 CONTENTS ix PART 2. STAGE TRANSITIONS AND DEMOGRAPHICS 183 Chapter 12. The Transition Through a Stage 185 12.1 12.2 12.3 12.4 12.5 12.6 12.7 12.8 12.9 The Sojourn Function Mean Sojourn Time, Expected Exit Age, and Expectation of Life The Variance of the Sojourn Time, Moments and Central Moments Remaining Sojourn Time and Its Expectation Fixed Stage Durations Per Capita Exit Rates (Mortality Rates) Exponentially Distributed Stage Durations Log-Normally Distributed Stage Durations A Stochastic Interpretation of Stage Transition Bibliographic Remarks Chapter 13. Stage Dynamics with Given Input 13.1 13.2 13.3 13.4 13.5 13.6 13.7 Input and Stage-Age Density The Partial Differential Equation Formulation Stage Content and Average Stage Duration Average Stage Age Stage Exit Rates 13.5.1 The Fundamental Balance Equation of Stage Dynamics 13.5.2 Average Age at Stage Exit Stage Outputs Which Recruitment Curves Can Be Explained by Cannibalism of Newborns by Adults? Bibliographic Remarks Chapter 14. Demographics in an Unlimiting Constant Environment 14.1 14.2 14.3 The Renewal Equation Balanced Exponential Growth The Renewal Theorem: Approach to Balanced Exponential Growth 185 187 189 190 197 199 201 202 206 209 211 211 212 217 219 221 222 224 226 230 237 239 240 241 244 Chapter 15. Some Demographic Lessons from Balanced Exponential Growth 255 15.1 15.2 15.3 15.4 Inequalities and Estimates for the Malthusian Parameter Average Age and Average Age at Death in a Population at Balanced Exponential Growth. Average Per Capita Death Rate Ratio of Population Size and Birth Rate Consequences of an Abrupt Shift in Maternity: Momentum of Population Growth Bibliographic Remarks Chapter 16. Some Nonlinear Demographics 16.1 16.2 A Demographic Model with a Juvenile and an Adult Stage A Differential Delay Equation Bibliographic Remarks 255 262 266 267 270 273 274 277 279 x CONTENTS PART 3. HOST–PARASITE POPULATION GROWTH: EPIDEMIOLOGY OF INFECTIOUS DISEASES 281 Chapter 17. Background 283 17.1 17.2 Impact of Infectious Diseases in Past and Present Time Epidemiological Terms and Principles Bibliographic Remarks Chapter 18. The Simplified Kermack–McKendrick Epidemic Model 18.1 18.2 18.3 A Model with Mass-Action Incidence Phase-Plane Analysis of the Model Equations. The Epidemic Threshold Theorem The Final Size of the Epidemic. Alternative Formulation of the Threshold Theorem Chapter 19. Generalization of the Mass-Action Law of Infection 19.1 19.2 19.3 Population-Size Dependent Contact Rates Model Modification The Generalized Epidemic Threshold Theorem Chapter 20. The Kermack–McKendrick Epidemic Model with Variable Infectivity 20.1 A Stage-Age Structured Model 20.2 Reduction to a Scalar Integral Equation Bibliographic Remarks Chapter 21. SEIR (→ S) Type Endemic Models for "Childhood Diseases" 21.1 21.2 21.3 21.4 The Model and Its Well-Posedness Equilibrium States and the Basic Replacement Ratio The Disease Dynamics in the Vicinities of the Disease-Free and the Endemic Equilibrium: Local Stability and the Interepidemic Period Some Global Results: Extinction, Persistence of the Disease; Conditions for Attraction to the Endemic Equilibrium Bibliographic Remarks Chapter 22. Age-Structured Models for Endemic Diseases and Optimal Vaccination Strategies 22.1 A Model with Chronological Age-Structure 22.2 Disease-Free and Endemic Equilibrium: the Replacement Ratio 22.3 The Net Replacement Ratio, and Disease Extinction and Persistence 22.4 Cost of Vaccinations and Optimal Age Schedules 284 289 291 293 293 295 297 305 305 306 307 311 311 313 316 317 318 321 325 332 339 341 341 348 351 358 CONTENTS 22.5 xi Estimating the Net Replacement Ratio: Average Duration of Susceptibility and Average Age at Infection. Optimal Vaccination Schedules Revisited Bibliographic Remarks 366 381 Chapter 23. Endemic Models with Multiple Groups or Populations 383 23.1 23.2 23.3 23.4 23.5 23.6 23.7 The Model Equilibrium Solutions Local Asymptotic Stability of Strongly Endemic Equilibria Extinction or Persistence of the Disease? The Basic Replacement Matrix, Alias Next-Generation Matrix The Basic Replacement Ratio as Spectral Radius of the Basic Replacement Matrix Some Special Cases of Mixing Bibliographic Remarks 384 388 394 399 404 406 411 416 PART 4. TOOLBOX 419 Appendix A Ordinary Differential Equations 421 A.1 A.2 A.3 A.4 A.5 Conservation of Positivity and Boundedness Planar Ordinary Differential Equation Systems The Method of Fluctuations Behavior in the Vicinity of an Equilibrium Elements of Persistence Theory Bibliographic Remarks Global Stability of a Compact Minimal Set Hopf Bifurcation Perron–Frobenius Theory of Positive Matrices and Associated Linear Dynamical Systems Bibliographic Remarks 421 424 428 433 436 441 442 444 Appendix B Integration, Integral Equations, and Some Convex Analysis 453 A.6 A.7 A.8 B.1 B.2 B.3 B.4 B.5 B.6 B.7 The Stieltjes Integral of Regulated Functions Some Elements from Measure Theory Some Elements from Convex Analysis Lebesgue–Stieltjes Integration Jensen's Inequality and Related Material Volterra Integral Equations Critical and Regular Values of a Function Bibliographic Remarks Appendix C Some MAPLE Worksheets with Comments for Part 1 C.1 C.2 Fitting the Growth of the World Population (Figure 3.1) Periodic Modulation of Exponential Growth in Closed Populations (Figures 3.2 and 3.3) 446 451 453 465 472 475 483 486 490 491 493 493 496 xii CONTENTS C.3 C.4 C.5 C.6 Fitting Sigmoid Population-Growth Curves (Figures 6.1 and 6.2) Fitting Bernoulli's Equation to Population Data of Sweden (Figure 6.3) Illustrating the Allee Effect (Figures 7.2–7.4) Dynamics of an Aquatic Population Interacting with a Polluted Environment (Figure 10.3E) 498 507 510 513 References 519 Index 537 Preface Now we only see models, like reflections in a mirror; but then we shall see face to face. Now I only know partially; but then I shall know as fully as I am myself known. St Paul, 1st letter to the Corinthians, 13:12 Plato's allegory of the cave suggests that we only see the shadows of reality, St Paul declares our knowledge patchwork, and Immanuel Kant distinguishes between the phenomenon we can recognize and the underlying noumenon (Ding an sich), which remains concealed. In more modern terms, we only see reality through models, first through the spontaneous models created by our senses, then through the deliberate models of science and arts, and finally through the meta-models of philosophy and theology. In this scenario, mathematics has its place as a form of symbolic modeling, namely, representation and analysis of reality through mathematical symbols and concepts. Biology, the science of life, has developed its own (nonmathematical) models, but recently the formulation of the dynamics of biological populations in mathematical equations, the analysis of these equations, and the reinterpretation of the results in biological terms has become a valuable source of insight, a selected sample of which has been collected in this book. As with any form of human cognition, mathematical modeling is imperfect, and its strengths and weaknesses must always be accounted for; they are discussed in more detail at the outset of our exposition. A treatise in mathematical biology can be organized according to biological topics or mathematical concepts and tools. We have chosen the first organizing principle; we start in Part 1 with unstructured single-species population models and then add the most rudimentary stage structure with variable stage duration. In Part 2, we develop the theme of stage transition more fully and apply it to a particular stage, human life, leading to the age-structured models of demographics. In Part 3, we consider the dynamic interplay of host and parasite populations, i.e., epidemic and endemic models. The theme of stage structure recurs in the form of the different stages of infection, and age-structure again plays a role in optimizing vaccination strategies. While the details can easily be seen from the table of contents, let us highlight some special features of this book. xiv PREFACE Realistic models should take account of seasonal and stochastic effects, but the incorporation of either of them makes it much more difficult to obtain mathematical results with a succinct biological interpretation. For that reason they are not dealt with in many mathematical biology books. Stochastic models (as well as statistical methods) are not included in this book either, but seasonality is at least incorporated in our treatment of single-species population dynamics. We develop concepts like birth-date dependent life expectancies and emigration probabilities. For an environment that does not limit population growth, we study the impact of time-periodic immigration on the development of a population the size of which would otherwise decrease. We also touch on the topic of seasonality for density-dependent growth, introducing the concept of asymptotic equality for population sizes. A great deal of care is spent in Part 1 of this book in deriving the classical single-species models linked to the names of von Bertalanffy, Verhulst, Beverton–Holt, Smith, Ricker, Gompertz, and Allee in the framework of continuous-time models. They are of fundamental importance because they often serve as building blocks in complex models, and knowing their origins is paramount for their correct use. This point is elaborated on when exploring the interaction of an aquatic population with a toxicant, where revisiting the Beverton–Holt growth function allows us to specify the mechanisms by which the toxicant affects individual growth. I prefer to discuss this planar model (and the beautiful and powerful phase-plane methods that come with it) in place of the traditional predator–prey and competition models; excellent treatments of the latter have been and are still being written by people who are more expert in this topic than me.As a second paradigm, a single-species model with two stages is investigated. There is a considerable amount of literature on stage-structured models with variable stage length, but apparently the most elementary model consisting of two ordinary differential equations has been overlooked. Surprisingly, it comes with bifurcation phenomena as complex as they can get in two space dimensions. The two planar models are qualitatively different, as the first can display multiple extinction equilibria (or, mathematically speaking, boundary steady states), while the second can have multiple persistence equilibria (interior steady states). The second part of the book fully develops the theme of stage structure in an age-dependent context, and demographic concepts such as life expectancy, expectation of remaining life, variance of life length, and their dynamic consequences get a mathematically rigorous treatment. They will later be linked to age-dependent models for infectious childhood diseases. Reflecting my personal research history, this book contains a relatively large part on host–parasite dynamics, i.e., epidemics and endemics of infectious diseases. A section on the Kermack–McKendrick epidemic model derives estimates from below and above for the final size of an epidemic that cannot be found elsewhere; they are sharp for basic replacement ratios that are either close PREFACE xv to 1 or very large. For childhood diseases we derive the formula for the length of the interepidemic period. The relation between the basic replacement ratio and the average age at infection (or, more correctly, average duration of susceptibility) is established without any restrictions on the survivorship. Sexual-disease transmission has highlighted the importance of population heterogeneity and the importance of multi-group models; in this book the relation between overall basic replacement ratio on the one hand and group-to-group basic replacement ratios on the other hand is interpreted in the light of the Perron–Frobenius theory for positive matrices. The interaction between population structure and dynamics has always fascinated me. This book concentrates on age-structure, stage structure and multiple diverse groups; the theory of general physiological structure is still under development and is not included. This text has grown from a series of courses I taught at Arizona State University over the years, ranging from the junior undergraduate to the advanced graduate level. So, hidden underneath the biological organization is an arrangement of the material according to mathematical difficulty, which the table of contents may not reveal clearly. The first four chapters of Part 1 mostly require one-variable calculus. A course in elementary differential equations may be helpful, but the methods of integrating factors, variation of parameters, separation of variables, and transformation of variables are developed from scratch and applied to more examples and practiced in more exercises than in a typical textbook for elementary differential equations. In Chapters 7–9, some advanced calculus may be needed; Chapter 9, on difference equations, is actually a celebration of the intermediate value theorem. Chapters 10 and 11 require an intermediate or advanced course in ordinary differential equations, notably phaseplane analysis and Hopf bifurcation in the plane. The required theorems are collected in a "toolbox." Part 2, on stage transition, again starts at the one-variable calculus level and provides ample material to revisit the definition of the integral, guided by a biological interpretation, and its properties. In a few instances, the concept of Stieltjes integration is needed; a full exposition of the Stieltjes integral of regulated functions is presented in the toolbox. In the last section, advanced calculus is again helpful. The epidemics division of Part 3 (Chapters 17–20) returns to the calculus level; Chapter 21 on childhood diseases requires an intermediate course on ordinary differential equations. Chapter 22, which approaches the question of optimal age-dependent vaccination strategies via linear optimization, has the most demanding background material, some measure theory and Lebesgue integration and elements of convex analysis which are also included in the toolbox. The background material for Chapter 23 (multi-group endemic models) is less demanding, but also less likely to be taught in a graduate course, namely the Perron–Frobenius theory of quasi-positive matrices and their role in systems xvi PREFACE of ordinary differential equations, the importance of which cannot be overestimated for models with an extended population structure. The more advanced last sections of both Part 1 and Part 3 are permeated by persistence theory, a dynamical systems method that allows us to determine which parts of an ecological or epidemic system do not become extinct on the one hand and remain bounded on the other hand. In spite of its abstract nature, the concept of a semiflow is used because it is the most natural in this context. Apart from this, the presentation in the toolbox is kept as elementary as possible by using uniform weak persistence as an intermediate step. Uniform weak persistence can often be checked by easy ad hoc considerations and implies the stronger property of uniform strong persistence under natural extra conditions. Another recurrent technique is the use of time-scale methods, which come in the form of quasi-steady-state simplifications, power expansions, and other approximations. They allow us to arrive at simpler models after starting from more complex ones which are formulated from first principles. They also help us to obtain simple approximation formulas which have a suggestive biological interpretation, in the context of not only differential equations but also of integral equations. I feel that scaling methods can best be learned by ad hoc applications, the paradigms in this being: the derivation of the Beverton–Holt, Ricker, and Allee population-growth functions; the discussion of the kind of recruitment functions that can be explained by cannibalism; the derivation of the formulas for the interval between outbreaks in childhood diseases; the relation between basic replacement ratio and average duration of the susceptible stage; and the complex formation approach to modeling infectious contacts in a multi-group population. Each chapter has exercises at its end, a few of them with solutions. Bibliographic remarks which provide further background information or make suggestions for further reading are added at the end of a chapter or after several shorter chapters. Chapter C in the toolbox consists of a collection of commented Maple worksheets which I used for producing some of the figures in the text. Acknowledgments This text grew from courses I taught at Arizona State University (ASU) in six semesters from Fall 1989 to Spring 1999. I thank the students who attended these courses and contributed to the development of the material in various ways, among them Zhilan Feng, Tim Lant, Irakli Loladze, Jimmy Mopecha, Natalia Navarova, and Jinling Yang. The only material never taught in a course at ASU is that in Chapter 11; it evolved from one of the lectures I gave at the Third Winterschool on Population Dynamics in Woudschoten (the Netherlands) in January 2000, which was organized by Odo Diekmann, Hans Heesterbeek, Fleur Kelpin and Bob Kooi. PREFACE xvii This book could not have been accomplished without the support of my colleagues in the Department of Mathematics and Statistics, starting with Tom Trotter, chairman at that time, who allowed me to teach one of the earlier courses with only two students, so that the project could get off the ground. The members of our Mathematical Biology group, Steve Baer, Frank Hoppensteadt, Yang Kuang, Hal Smith, and the late Betty Tang, created the nourishing and inspiring environment in which such an endeavor could flourish. Steve Baer spent many hours with me in performing the AUTO calculations for Figures 11.5 and 11.6 in Chapter 11, and Hal Smith read parts of the typescript. John McDonald helped with the convex optimization material in Toolbox B.3, Matt Hassett and Dennis Young with the probability distributions in Chapter 13, and Joe Rody was always willing to answer my Maple questions. Renate Mittelmann and Jialong He and their student assistants tended to the various computers which were used over the years in producing this book. I would like to mention a few things that influenced this book without leaving direct evidence in the contents: the papers by Frank Hoppensteadt and Paul Waltman (1970, 1971) on an epidemic model suggested by Ken Cooke (1967) with which my "Doktorvater", Willi Jäger, made me acquainted as preprints and which got me started as a mathematical biologist; the book by Krasnosel'skii, Positive Solutions of Operator Equations (1964), the methods of which made my Ph.D. dissertation fly (they actually play a role in Chapter 23); the discussions on mathematical biology in Willi Jäger's research group, particularly with Wolfgang Alt and Konrad Schumacher; the joint work with Odo Diekmann, Mats Gyllenberg, Henk Heijmans, and Hans Metz on physiologically structured population models stretching over two decades; and, at ASU, the collaboration with Hal Smith. I am grateful for the many detailed and constructive reviews the typescript has obtained, from Reinhard Bürger, Odo Diekmann, Herb Hethcote, and anonymous reviewers. Many others have helped by commenting on the typescript, providing important references, or through their general support, among them Fred Brauer, Carlos Castillo-Chavéz, John Chance, Jim Collins, Jim Cushing, Patrick DeLeenheer, Klaus Dietz, Zhilan Feng, Bob Kooi, Bas Kooijman, Karl Hadeler, Mimmo Iannelli, Ulrich Krause, Torsten Lindström, Maia Martcheva, Fabio Milner, Roger Nisbet, Diana Thomas, and Pauline van den Driessche. As part of their honors projects, Michael Belisle, Carrie Eggen, and Sumir Mathur read some sections and worked some exercises and tested some of the Maple worksheets for snags. During most of the years I worked on this book, I had continuing partial support from the Division of Mathematical Sciences at the National Science Foundation, under the program directorate of Mike Steuerwalt. It must have been in the early 1990s that I first mentioned my germinating book project to Simon Levin and showed him a few sections of an early version. He immediately became supportive, with a view to later publication. I warned xviii PREFACE him, though, that it might take quite a few more years (it turned out to be more than 10) and asked him when he would lose interest. He answered "probably never," and he never did, eventually guiding the project to Princeton University Press. Thanks also to David Ireland and his colleagues at Princeton University Press, and Sam Clark, copy-editor and typesetter at T&T Productions Ltd, who lent their experience and hard work to the task of bringing the typescript to hard- and soft-cover reality. I dedicate this book to my wife Adelheid, for her love and her support. During the years in which this book took shape, she kept the family and the house running (along with her own professional career), and shielded me from the time-consuming and annoying disturbances of everyday life. She also checked the typescript for nonmathematical mistakes and typos; all remaining errors are my responsibility, though. I thank my daughters, Ruth and Clara, for their tolerance. Last, above all, praise to the Holy Spirit, source of knowledge and persistence, for granting the completion of this book. Mathematics in Population Biology Chapter One Some General Remarks on Mathematical Modeling ` κ µ´ ρoυς γ ὰρ γ ιν ώσ κoµν St Paul, 1st letter to the Corinthians, 13:9 Modeling is an attempt to see the wood for the trees. A model is a simplification or abstraction of reality (whatever that is), separating the important from the irrelevant. Actually, modeling is a part of our existence. If we want to be philosophical, we could say that we do not perceive reality as it is, but only realize a model our mind has designed from sensory stimuli and their interpretation. It seems that certain animal species perceive different models of reality which, compared to ours, are based more on hearing and smell than on sight. Many philosophers have had much deeper thoughts on this problem than I present here; following Plato's famous allegory of the cave we may say that we only see the shadows of reality, or, following Kant, that we see the phenomena rather than the noumena. St Paul (1st letter to the Corinthians, 13:9), puts it this way: We obtain our knowledge in parts, and we prophesy in parts. St Paul, of course, had a broader and deeper reality in mind than we are concerned with here, but if we apply his statement to a more restricted reality, we may (rather freely) paraphrase it as follows: We obtain our knowledge from models, and we make our predictions on the basis of models. Since, in a wide sense, we are modeling anyway, modeling in the strict sense is the purposeful attempt to replace one model (the so-called "real world," which we typically accept without questioning) by another, deliberate, model which may give us more insight. Essentially, there seem to be two incentives for doing so: either the "real-world" model is too complex to obtain the desired insight and so is replaced by a simpler or more abstract one. Or the "real-world" model does not allow certain experiments for ethical, practical or other reasons and is replaced by a model in which all kinds of changes can be readily made and their consequences studied without causing harm. The word model presumably traces back to the Latin word modulus, which means "little measure" (Merriam-Webster, 1994), alluding to a small-scale physical representation of a large object (e.g., a model airplane). 2 CHAPTER 1 Mathematical modeling uses symbolic rather than physical representations, unleashing the power of mathematical analysis to increase scientific understanding. It can conveniently be divided into three stages (cf. Lin, Segel, 1974, 1988). (i) Model formulation: the translation of the scientific problem into mathematical terms. (ii) Model analysis: the mathematical solution of the model thus created. (iii) Model interpretation and verification: the interpretation of the solution and its empirical verification in terms of the original problem. All three steps are important and useful. Already the first step—model formulation—can lead to considerable insight. For building a mathematical model, one needs clear-cut assumptions about the operating mechanisms, and it often turns out to be an unpleasant surprise that the old model—the "real world"—is far less understood than one thought. In many cases the modeling procedure—at least if one chooses parameters that are meaningful—already teaches what further knowledge is needed in order to apply the mathematical model successfully; the model analysis and its interpretation help to determine to what extent and precision new information and new data have to be collected. Analytic and numerical tools allow the extrapolation of present states of the mathematical model into the future and, sometimes, into the past. Assumptions, initial states, and parameters can easily be changed and the different outcomes compared. So, models can be used to identify trends or to estimate uncertainties in forecasts. Different detection, prevention and control strategies can be tested and evaluated. While the model analysis may require sophisticated analytic or computational methods, mathematical modeling ideally leads to conceptual insight, which can be expressed without elaborate mathematics. A model is, purposefully, a simplification or abstraction; very often it is an oversimplification or overabstraction. Insight obtained from a model should be checked against empirical evidence and common sense. It can also be checked against insight from other models: how much does the model's behavior depend on the degree of complexity, on the form of the model equations, on the choice of the parameters? Dealing with a concrete problem, a modeler should work with a whole scale of models starting from one which is as simple as possible to obtain some (hopefully) basic insight and then adding complexity and checking whether the insight is confirmed. Starting with a complex model has the risk that the point will be missed because it is obscured by all the details. The use of a whole range of models also educates the modeler on how critically qualitative and quantitative results depend on the assumptions one has made. A mathematical modeler should become a guardian keeping him/herself and others from jumping to premature conclusions. SOME GENERAL REMARKS ON MATHEMATICAL MODELING 3 When modeling concrete phenomena, there is typically a dilemma between incorporating enough complexity (or realism) on the one hand and keeping the model tractable on the other. In some cases the modeler feels forced to add a lot of detail to her/his model in order to grasp all the important features and ends up with a huge number of unknown parameters, which may never be known even approximately and cannot be determined by fitting procedures either. Such a mathematical model will be of limited value for quantitative and maybe even qualitative forecasts, but still has the other benefits described above. Mathematical modeling has its place in all sciences. This book, according to the interests of the author, concerns mathematical models in the biosciences, or rather a very small area within the biosciences, namely population biology. Another restriction is the one to deterministic models (as opposed to stochastic models), which neglect the influence of random events. To some degree one can dispute whether stochastic models are more realistic than deterministic models; there is still the possibility that everything is deterministic, but just incredibly complex. In this case, stochasticity would simply be a certain way to deal with the fact that there are many factors we do not know. Stochastic models are theoretically more satisfying as they count individuals with integer numbers, while deterministic models, usually differential or difference equations, have to allow population sizes that are not integers. (This actually is of no concern if population size is modeled as biomass and not as number of individuals, but the latter is not always appropriate.) While a typical tool of deterministicmodel analysis consists of discussing large-time limits, stochastic models take account of the truism that nothing lasts forever and make it possible to analyze the expected time until extinction—a concept that has no counterpart in deterministic models. In many cases, deterministic models can theoretically be justified as approximations of stochastic models for large populations sizes; however, the population size needed to make the approximation good enough may be unrealistically large. Nevertheless, deterministic models have the values which we described above, as long as one keeps their limitations in mind. The latter particularly concerns predictions, which are of very limited use in this uncertain world if no confidence intervals for the predicted phenomena are provided. Bibliographic Remarks Other discussions of mathematical modeling can be found in Lin, Segel (1974, 1988), Hethcote,VanArk (1992), Kooijman (2000), and Kirkilionis et al. (2002). Section 1.1 in Lin and Segel's On the Nature of Applied Mathematics more generally speaks about the scope, purpose, and practice of applied mathematics (also in comparison with pure mathematics and theoretical science), while Section 1.6 in Hethcote, Van Ark (1992), Purposes and limitations of epidemiological modeling, refers to the modeling of infectious diseases and of the 4 CHAPTER 1 AIDS epidemic in particular. Their respective remarks, however, easily specialize or generalize to mathematical modeling. A "playful" approach to the role of modeling and mathematics in the sciences can be found in Sigmund (1993, Chapters 1, 3, and 9). The idea that a concrete problem should be addressed by a chain of models has been realized in a textbook by Mesterton-Gibbons (1992), who calls this the layered approach. Section 1.2 in Kooijman (2000) shares the philosophical touch of this chapter; the outlooks are rather similar for modeling itself, but quite different in the meta-modeling aspects. For a discussion of deterministic versus stochastic models see the papers by Nåsell (n.d.) and Section 1.6 in the book by Gurney, Nisbet (1998). Genetic and evolutionary aspects of population biology are not touched on at all in this book, and I refer to the books by Bürger (2000), Charlesworth (1980, 1994), Farkas (2001), and Hofbauer, Sigmund (1988, 1998) as sources for their mathematical modeling and analysis. PART 1 Basic Population Growth Models Chapter Two Birth, Death, and Migration In a natural environment, the dynamics of one species can hardly be isolated from those of other species. The particular species we are interested in may depend on others as resources or may serve as a resource itself (prey–predator relation), and may compete or cooperate with other species in the quest for food and habitat. Nevertheless, dividing to conquer, mathematical modeling first concentrates on one population and expresses the change of its size in terms of concepts like birth rates, mortality rates, and emigration and immigration rates which take the form of either population rates or per capita rates. 2.1 The Fundamental Balance Equation of Population Dynamics As a first goal, mathematical models try to describe the size of a population as a function of time. A second goal may be to describe the distribution of the population over space or according to important individual characteristics like age, body size, protein content, etc., but we ignore population structure at this point. So, let N (t) denote the size of the single-species population we are focusing on at time t. Depending on the species, the population size may be the number of individuals or the collective biomass of the population. It can also be the average spatial density of the population or, for aquatic populations, the average concentration per water unit. The rate of change of the population size is given by the time derivative N (t) = dN (t)/dt. If we consider a species of mammal, we may like to have the population size to be given in numbers of individuals. Then the size can change by four processes: births, deaths, emigration, and immigration. If there is no emigration or immigration, we speak of a closed population, otherwise of an open population. The fundamental balance equation of population dynamics takes these four processes into account, N (t) = B(t) − D(t) + I (t) − E(t), (2.1) where, at time t, B(t) denotes the population birth rate, D(t) the population death rate, I (t) the immigration rate, and E(t) the emigration rate. More precisely, B(t) is the number of births per time unit, at time t, D(t) is the number of deaths per unit of time, etc. This balance equation can equivalently be written 8 CHAPTER 2 as t N (t)−N (r) = t B(s) ds − r t D(s) ds + r t I (s) ds − r E(s) ds, t > r. r As an illustration, we give the example of Sweden in 1988 (Preston et al., 2001, Box 1.1, p. 3). There t is January 1, 1989, and r is January 1, 1988, and N (t) = 8 461 554, N (r) = 8 416 599. The number of births and deaths in 1988 are t t B(s) ds = 112 080, D(s) ds = 96 756, r r while the number of immigrants and emigrants are t t I (s) ds = 51 092, E(s) ds = 21 461. r r If the individuals are very small, it may be better to describe the population size in terms of biomass. Then B would be the biomass acquisition rate through food uptake and D the biomass loss rate due not only to death but also to nonfatal starvation. As for births, deaths, and emigration, it also makes sense to speak about per capita rates. The per capita birth rate, β(t), is the average number of offspring born to one typical individual per unit of time, at time t. At this point, however, the per capita birth rate is just a theoretical construct. In sexually reproducing populations, for example, only females give birth, and this only in a certain agewindow, and only after mating. So, in order to determine the per capita birth rate, one would need to know the sex ratio, the age-distribution of the female population, and the mating probability in order to determine the per capita birth rate as defined above. To complicate matters even more, all these ingredients may change with time. As for the per capita death or mortality rate, µ(t), even the definition takes substantial effort. Let Π (t, r), t r, denote the probability of surviving from time r to time t, Π (r, r) = 1, 0 Π (t, r) 1. The probability of dying from time r to time t is 1 − Π (t, r). The per capita mortality rate µ(t) is defined as µ(t) = lim h→0 1 − Π (t + h, t) , h assuming that these limits exist. The number of people dying from time t to time t + h is N (t)[1 − Π (t + h, t)], so the population death rate at time t is N (t)[1 − Π (t + h, t)] = µ(t)N (t). h→0 h D(t) = lim 9 BIRTH, DEATH, AND MIGRATION Π can be recovered from µ if the probability of surviving from r to t is independent from the past up to time r, in particular if mortality does not depend on age. Then the following multiplicative rule holds: Π (t, r) = Π (t, s)Π (s, r), t s r, Π (r, r) = 1. This implies Π (t + h, r) − Π (t, r) 1 − Π (t + h, t) =− Π (t, r), h h and the following differential equation holds for Π , d Π (t, r) = −µ(t)Π (t, r), dt Π (r, r) = 1. We fix r and set p(t) = Π (t, r). In terms of p, the differential equation can be solved by separation of variables, −µ(t) = p (t) d = ln p(t), p(t) dt p(r) = 1. We integrate both sides and use the fundamental theorem of calculus, t µ(s) ds = ln p(t) − ln p(r) = ln p(t) − ln 1. − r Exponentiating both sides provides t µ(s) ds . Π (t, r) = p(t) = exp − (2.2) r 2.2 Birth Date Dependent Life Expectancies For later use, let us calculate the life expectancy, L(r), of an individual born at time r. Then r + L(r) is the expected time of death. We first assume that there is a maximum finite life span strictly less than some number c. This means that the individual born at time r will be dead by time r + c, i.e., p(r + c) = Π (r + c, r) = 0. We choose a partition r = t1 < · · · < tn+1 = c + r of the interval [r, r + c] such that tj +1 − tj is very small for all j . Observe that pj = p(tj +1 ) − p(tj ) = Π (tj +1 , r) − Π (tj , r) is the probability that the individual dies between tj and tj +1 . In this case, the time of death is some number sj between tj and tj +1 . The expected time of death is given as the limit of Riemann–Stieltjes sums, L(r) + r = lim n j =1 sj pj = lim n j =1 sj [p(tj ) − p(tj +1 )], 10 CHAPTER 2 where the lengths of all intervals [tj , tj +1 ] tend to 0. Taking the sums apart and changing the index of summation we obtain L(r) + r = lim n sj p(tj ) − j =1 n+1 sj −1 p(tj ) . j =2 We set s0 = r and sn+1 = r + c. Then r = s0 < · · · < sn+1 = r + c is a partition of r + c and L(r) + r can be expressed as a limit of Riemann sums, i.e., as a Riemann integral, L(r) + r = lim n+1 (sj − sj −1 )p(tj ) − sn+1 p(tn+1 ) + s0 p(t1 ) j =1 r+c p(t) dt − (r + c)p(r + c) + rp(r). = r Since p(r) = 1 and p(r + c) = 0, by substituting t = s + r, c c L(r) = p(s + r) ds = Π (s + r, r) ds. 0 0 If there is no finite maximum life span, we let c → ∞, and L(r) is the improper Riemann integral ∞ L(r) = Π (s + r, r) ds. 0 In terms of per capita mortality rates µ, by (2.2), r+s ∞ ∞ L(r) = exp − µ(t) dt ds = exp − 0 0 r s µ(t + r) dt ds. 0 In the important special case that the per capita mortality rate is constant, L(r) does not depend on r and the life expectancy is given by ∞ 1 L= (2.3) e−µs ds = . µ 0 Mathematically, there is not much difference between leaving the population by death or by emigration. The per capita emigration rates are denoted by η(t) and can be defined from the probabilities Ξ (t, r) of not emigrating from time r to time t. A derivation similar to the one above gives that the mean sojourn time of an individual, born at time r, in this particular population (death neglected) is given by ∞ 0 s exp − η(t + r) dt ds, 0 and by 1/η if the per capita emigration rate is constant. If both death and emigration are regarded and the per capita mortality and emigration rates are 11 BIRTH, DEATH, AND MIGRATION constant, the mean sojourn time in the population is then 1/(µ + η), because µ + η is the rate of leaving the population by either death or emigration. 2.3 The Probability of Lifetime Emigration Again for later use, let us determine the probability that an individual born at time r will emigrate during its lifetime. Fix the birth date r and let ξ(t) = Ξ (t, r) denote the probability of not having emigrated until time t. The probability of emigrating in the time interval between tj −1 and tj , where tj > tj −1 r, is ξ(tj −1 ) − ξ(tj ). So, the unconditional probability of emigrating between those times is [ξ(tj −1 ) − ξ(tj )]p(sj ) with some sj between tj −1 and tj . Here p(s) again denotes the probability to be still alive at time s. We assume that ξ is continuously differentiable as is the case when we have a continuous emigration rate η = ξ . By the mean value theorem, there exists some s̃j between tj −1 and tj such that the probability of emigrating between these times is −ξ (s̃j )(tj − tj −1 )p(sj ). The probability of emigrating before some time T > r is obtained by choosing a partition r = t0 < · · · < tn = T , by summing these probabilities, and letting the number of partition points tend to infinity: lim n→∞ n (−1)ξ (s̃j )p(sj )(tj − tj −1 ) = − j =1 T ξ (t)p(t) dt. r The probability of emigrating during one's lifetime, Pout (r), is obtained by taking the limit T → ∞. If we have continuous per capita mortality and emigration rates, µ(·) and η(·), then t t p(t) = exp − µ(s) ds and ξ(t) = exp − η(s) ds r r by (2.2). So, the probability of emigrating during one's lifetime is t ∞ η(t) exp − [η(s) + µ(s)] ds dt. Pout (r) = r r If η and µ are constant, so is this probability, and Pout = η η+µ (2.4) is the product of the mean sojourn time in the population and the constant per capita emigration rate. The following relations hold between population and per capita rates: B(t) = β(t)N (t), D(t) = µ(t)N (t), E(t) = η(t)N (t). 12 CHAPTER 2 It makes little sense in this context to speak about a per capita immigration rate; for it may well happen that N (t) = 0 and I (t) > 0, if a species invades a formerly empty habitat. We mention that these considerations also apply if we only consider the part of a population in a certain stage (e.g., larvae). The birth rate has then to be interpreted as a stage influx rate (stage input, the hatching rate for larvae, for example) and the emigration rate as the stage exit rate (stage output, the pupation rate for larvae, for example). For constant per capita rates, 1/η gives the average stage length (deaths disregarded), and 1/(µ + η) the stage sojourn time, while η/(η + µ) is the probability of surviving the stage and µ/(η + µ) the probability of dying during the stage (cf. Chapter 11). Chapter Three Unconstrained Population Growth for Single Species If we look at a single species exclusively, the population rates in the fundamental balance equation of population dynamics would only depend on the size of the respective population, N , and on the time variable. Although in most cases the dynamics of this population cannot be isolated from the dynamics of others, it is still useful to look at single-species dynamics for a number of reasons: • Single-species dynamics occur in laboratory environments. • Sometimes it may be possible to blend other species into the environment. • Understanding one-species dynamics helps us to understand the more complex multi-species dynamics. There are two separate cases to consider: unconstrained population growth, where the per capita rates depend only on environmental conditions, i.e., they are functions of time only or are just constant; and density-dependent population growth, where the per capita rates are affected by the population size, N . The first case occurs if a species has unlimited resources and is not subject to predators or competitors; the per capita rates are then just functions of time t or may even be constant. 3.1 Closed Populations If there is neither immigration nor emigration, equation (2.1) takes the form of a linear ordinary differential equation, N = ρ(t)N, ρ(t) = β(t) − µ(t), (3.1) N (t0 ) = N0 , where N0 is the population size at some time t0 (e.g., the time when we start observing the population). N0 and t0 are called the initial data of our model. We assume that the per capita birth rate β and per capita mortality rate µ are continuous. Their difference ρ is sometimes called the intrinsic (per capita) growth rate. 14 CHAPTER 3 Equation (3.1) can be solved by a method called separation of variables, where we collect all terms containing the dependent variable N on one side of the equation: N (t) d ρ(t) = = ln N (t). N (t) dt We integrate both sides and use the fundamental theorem of calculus, t t t0 ln N (t) − ln N (t0 ) = ρ(s) ds = ρ(s) ds − ρ(s) ds. 0 t0 0 Exponentiating both sides provides Q(t) , Q(t0 ) t Q(t) = exp ρ(s) ds . N (t) = N0 (3.2) 0 Q is called a fundamental solution of (3.1). Notice that, in the formula for Q, the lower integration limit 0 can be replaced by a fixed, but arbitrary, number r0 , which can be chosen at convenience. Actually, Q can be chosen as the exponential of any anti-derivative of ρ. If the per capita birth and death rates are constant, so is ρ(t) = ρ̄, and we obtain N (t) = N (t0 )eρ̄(t−t0 ) . This is the notorious exponential growth, in which a population doubles (or halves) its size in a fixed amount of time (see Exercise 3.1) if ρ̄ > 0 (or ρ̄ < 0, respectively). The British economist Thomas Malthus (1798) was one of the first to be concerned about its consequences, and the parameter ρ̄ is often called the Malthusian parameter in his honor. Other names for ρ̄ are intrinsic rate of natural increase or intrinsic growth constant. In Figure 3.1 we have fitted a line and an exponential curve to the world population data from 1950 to 1985 taken from Keyfitz, Flieger (1990, p. 105). The doubling time for the exponential curve is 36.34 years. For details see the Maple worksheet in Toolbox C. 3.1.1 The Average Intrinsic Growth Rate for Periodic Environments Constant per capita birth and mortality rates are highly unlikely for most natural populations; rather they are usually subject to seasonal fluctuations. A more realistic, though still rather idealistic, assumption is a periodic time-dependence where β(t) and µ(t), and so ρ(t), are periodic functions of t (with the same period), i.e., ρ(t) = ρ(t + p), t ∈ R. 15 UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES exponential fit linear fit 5 5 4 4 3 3 1950 1960 1970 1980 1990 1950 1960 1970 1980 1990 Figure 3.1. Linear and exponential fits to the growth of the world population (in billions) from 1950 to 1985 according to data from Keyfitz, Flieger (1990, p. 105). The fixed number p is the period of β, µ, and ρ. In many applications, p will be one year. Example. Let ρ(t) = 2 + sin t. We could find a solution of (3.1) by substituting ρ into (3.2). But let us instead do it from scratch. Separating the variables we obtain 2 + sin t = N d = ln N. N dt Integrating both sides we have ln N (t) = 2t − cos t + c with an integration constant c. Exponentiating, N (t) = ke2t e− cos t , k = ec . So, the population size exhibits exponential growth with rate ρ̄ = 2, modulated by the periodic function e− cos t . N , given in this form with an arbitrary constant k, is called a general solution. The constant k can be adjusted to make N satisfy the appropriate initial data. For instance, if N (1) = 5, then, by substitution into the general solution, 5 = ke2 e− cos 1 . Solving for k and fitting it into the general solution, N (t) = 5e2(t−1) ecos 1−cos t . The behavior of the population size in this example, periodically modulated exponential growth, makes us wonder whether this is a general phenomenon. More precisely, we ask ourselves the following question. If ρ is periodic, is there a number ρ̄ and a periodic function q (with same period p) such that (t−t0 )ρ̄ q(t) ? N (t) = N (t0 )e q(t0 ) 16 CHAPTER 3 By (3.2), N (t) = N (t0 )(Q(t)/Q(t0 )). We look for ρ̄ ∈ R such that q(t) = e−ρ̄t Q(t) is p-periodic. This is equivalent to σ (t) = ln q(t) being p-periodic. By the form of Q in (3.2), t σ (t) = ρ(r) dr − ρ̄t. 0 Obviously, σ (0) = 0. A minimum requirement for σ to be p-periodic is σ (p) = 0 as well. So, we set 1 p ρ̄ = ρ(r) dr. p 0 With this choice, σ inherits p-periodicity from ρ and so does q(t) = eσ (t) . Indeed, σ (t) = ρ(t) − ρ̄, so σ (t +p)−σ (t) = 0 for all t. Thus σ (t +p)−σ (t) is a constant function of t and σ (t + p) − σ (t) = σ (p) − σ (0) = 0 for all t. Notice that q(t) = eσ (t) is a bounded, strictly positive function. Let us summarize. For periodic ρ, we define ρ̄ to be the time average of ρ, 1 p ρ̄ = ρ(s) ds. (3.3) p 0 ρ̄ is the average intrinsic growth rate for periodic ρ. We have the following result. Theorem 3.1. Consider a closed population, N = ρ(t)N , where ρ is a pperiodic function and ρ̄ its time average. Then there exists some p-periodic function q such that the solution of N of (3.1) satisfies N (t) = N (t0 )e(t−t0 )ρ̄ q(t) . q(t0 ) Remark 3.2. In other words, the fundamental solution Q in (3.2) has the form Q(t) = eρ̄t q(t) with a p-periodic function q. q is continuously differentiable and strictly positive. With some justification one can again call ρ̄ the Malthusian parameter if ρ is periodic. Notice the following relation from Theorem 3.1, N (t0 + mp) = N (t0 )emρ̄ , m ∈ Z, which actually holds for all times t0 . This means that the population grows geometrically from season to season with growth factor eρ̄ . Theorem 3.1 is a special case of Floquet theory in the context of which ρ̄ is called a Floquet exponent and eρ̄ a Floquet multiplier. We call Q(t) = eρ̄t q(t) the Floquet representation of the fundamental solution Q. UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES 17 Illustration We consider a population subject to a periodic environment with a period of one year. We want to mimic the situation in which the per capita birth rate has a much more pronounced seasonality than the per capita mortality rate; in other words, the reproductive season is relatively short. Furthermore, we want the birth rate to be at its maximum when the mortality rate is at its minimum and vice versa, β̃(t) = (1 + cos(2πt))8 , µ̃(t) = 1 + 0.5 sin(2π t − π/2) (see Figure 3.2A). Both per capita rates have been normalized such that their averages are 1. We combine these normalized per capita birth and death rates into an intrinsic growth rate ρ in three different ways as follows. (A) ρA (t) = β̃(t) − µ̃(t), ρ̄A = 0. (B) ρB (t) = β̃(t) − 0.7µ̃(t), ρ̄B = 0.3. (C) ρC (t) = β̃(t) − 1.9µ̃(t), ρ̄C = −0.9. Figure 3.3 illustrates Theorem 3.1 for these various choices of ρ. If ρ̄ = 0, the population size remains in periodic motion; for ρ̄ = 0 the population increases or decreases exponentially with periodic modulation. The Maple worksheet used to produce Figures 3.2 and 3.3 can be found in Toolbox C. 3.1.2 The Average Intrinsic Growth Rate for Nonperiodic Environments If ρ is not necessarily periodic, we may consider the asymptotic time average 1 t→∞ t ρ̄ = lim t ρ(s) ds (3.4) 0 provided that this limit exists. If ρ is periodic, this definition is consistent with our previous one. If t > p, choose m ∈ N such that t ∈ [mp, (m + 1)p). Obviously, t → ∞ if and only if m → ∞. Let ρ̄ be defined by (3.3). We split up the integral, 1 t 0 t ρ(s) ds = m−1 1 (j +1)p 1 t ρ(s) ds + ρ(s) ds. t t mp jp j =0 18 CHAPTER 3 Per capita birth and death rates 5 A Immigration rate B birth 3 4 3 2 2 death 1 1 0 0 1 2 3 0 1 2 3 Figure 3.2. Fictitious periodic per capita birth, death, and immigration rates for a population with short reproductive season. The period is one year with the year beginning at the time when the birth rate is at its maximum. The immigration rate has its maximum at the same time, while the death rate is at its minimum. The birth and death rates have been normalized such that their averages are 1. Using substitution and the periodicity of ρ, 1 t 0 t m−1 1 p 1 t−mp ρ(s) ds = ρ(s + jp) ds + ρ(s + mp) ds t t 0 j =0 0 1 t−mp mp ρ̄ + ρ(s) ds. = t t 0 t−mp Since mp t < mp + p, we have (mp/t) → 1 and (1/t) 0 ρ(s) ds → 0 as t → ∞, which shows that ρ̄ in (3.4) is the same as in (3.3). So, we have proved the following result. Proposition 3.3. The asymptotic time average of a periodic function equals its time average. If the limit in formula (3.4) exists, we call ρ̄ the average intrinsic growth rate. Theorem 3.4. Assume that the average intrinsic growth rate ρ̄ in (3.4) is defined. Then the population increases exponentially if ρ̄ > 0, and decreases exponentially if ρ̄ < 0. Proof. Let ρ̄ > 0. Then there exists some r > 0 such that 1 t ρ(s) ds 21 ρ̄, t r. t 0 By (3.2), Q(t) eρ̄t/2 , t r, and N (t) increases exponentially. The case ρ̄ < 0 is left as an exercise. 19 UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES y(t) A B 1.4 1.2 1.0 0.8 0.6 0.4 0.2 0 20 15 10 5 0 2 4 t y(t) 0 6 8 0 10 2 4 t 6 8 10 C 1.4 1.2 1.0 0.8 0.6 0.4 0.2 0 0 2 4 t 6 8 10 Figure 3.3. Dynamics of closed populations under various combinations of the per capita birth and death rates of Figure 3.2. (A) ρ̄ = 0. (B) ρ̄ = 0.3. (C) ρ̄ = −0.9. 3.2 Open Populations If immigration and emigration are included, we have N (t) = ρ(t)N + I (t), ρ(t) = β(t) − µ(t) − η(t). (3.5) The differential equation in (3.5) is called a nonhomogeneous linear equation, while (3.1) is called a homogeneous linear equation. The equations in (3.5) can be solved by using an integrating factor. We collect the N terms on one side of the equation and multiply by some function F (t) (the integrating factor): F (t)N (t) − ρ(t)F (t)N (t) = F (t)I (t). If F (t) = −ρ(t)F (t), then the product rule allows us to write this equation as d (F (t)N (t)) = F (t)I (t). dt 20 CHAPTER 3 The first equation is of the form (3.1), so a special solution is t 1 . ρ(s) ds = F (t) = exp − Q(t) 0 Integrating the second equation gives F (t)N(t) − F (t0 )N (t0 ) = t F (s)I (s) ds. t0 We solve for N (t): t Q(t) Q(t) I (s) + ds, Q(t0 ) Q(s) t0 t Q(t) = exp ρ(r) dr . N (t) = N0 (3.6) 0 Again, in the formula for the fundamental solution Q, the lower integration limit 0 can be replaced by a fixed, but arbitrary, number r0 which can be chosen at convenience. In other words, we can choose Q as the exponential of any anti-derivative of ρ. Formula (3.6) is often called the variation of parameters formula or the variation of constants formula because it can also be found by a method called variation of parameters. The phenomenon that the nonhomogeneous equation can be solved in this way, in terms of the fundamental solution of the homogeneous equation and the nonhomogeneous term I , is also called Duhamel's principle. We notice from formula (3.6) that a solution to (3.5) is uniquely determined by its initial data N = N0 at t = t0 . Example. 1 1 N + 3, N (1) = N0 . t2 t Here the immigration rate tends to 0 as t → ∞, and the per capita death and emigration rates together are larger than the per capita birth rate, though the difference tends to 0 as well. As ρ(t) = −1/t 2 is not defined at t = 0, we do not use the formula for Q in (3.6), but take Q as the exponential of a convenient anti-derivative of ρ, Q(t) = e1/t . Then t 1 −1/s e ds. N (t) = N0 e(1/t)−1 + e1/t 3 1 s N = − Substituting s = −1/r, N (t) = N0 e(1/t)−1 + e1/t −1/t −1 (−1)rer ds. UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES 21 Finding an anti-derivative for rer (integrating by parts if we cannot guess it), r=−1/t N (t) = N0 e(1/t)−1 + e1/t [−rer + er ]r=−1 = N0 e(1/t)−1 + e1/t 1 −1/t e + e−1/t − e−1 − e−1 . t Simplifying, N (t) = N0 e(1/t)−1 + 1 + 1 − 2e(1/t)−1 . t We notice that N (t) → (N0 − 2)e−1 + 1 as t → ∞, i.e., the population size approaches a positive constant after a sufficiently long time. 3.2.1 Nonzero Average Intrinsic Growth Rate We now consider an open population in which the per capita birth rate is smaller than the per capita death and emigration rates together, in the asymptotic time average. It is in this situation that immigration has the biggest impact on the population development. An example is a fish pond in which the fish population cannot survive on its own and so is stocked regularly. We show that the solutions of (3.5) "forget" their initial data if the asymptotic time average ρ̄ in (3.4) is negative. Theorem 3.5. Assume that ρ̄ < 0 and consider two solutions N and M of (3.5) satisfying different initial data. Then N (t) − M(t) is a solution of (3.1) which converges to 0 as t → ∞. Proof. Let Ñ = N − M be the difference of the two solutions. Then Ñ = ρ(t)Ñ and Ñ (t) decreases exponentially to 0 by Theorem 3.4. We now return to the scenario in which the various per capita rates are periodic and assume that the immigration rate is also periodic with the same period p. Example. Let ρ(t) = cos(t) − 1, I (t) = 1 − cos(t). Rather than using the variation of parameters formula, we go through the integrating procedure again. The equation for the integrating factor F specializes to F = −(cos(t) − 1)F. Separating the variables we find the special solution F (t) = et−sin(t) . 22 CHAPTER 3 With this integrating factor the original differential equation becomes d d (F (t)N (t)) = F (t)I (t) = et−sin t (1 − cos(t)) = et−sin(t) . dt dt Integrating both sides, F (t)N (t) = et−sin(t) + c, and dividing by F (t), N (t) = 1 + cesin(t)−t . We see that N (t) → 1 as t → ∞. This asymptotic behavior is more regular than we can generally expect, as we will see in an example below. It is due to the fact that this example is of the form N = ρ(t)N − ρ(t)A with a positive constant A and that M(t) = A is a special solution. If the asymptotic time average ρ̄ exists and is strictly negative, it follows from Theorem 3.5 that N (t) − M(t) → 0 as t → ∞ and we have the following result. Notice that ρ does not need to be periodic. Corollary 3.6. Consider N = ρ(t)N − ρ(t)A with a constant A and assume that ρ̄ < 0. Then N (t) → A as t → ∞. In particular, if N = −αN + A with positive constants α and A, then N (t) → A/α as t → ∞. We will see a less regular asymptotic behavior in the next example. Example. Let ρ(t) = ρ̄ = −2 and again I (t) = 1 − cos t. We could specialize the general formula (3.6) or redo the integration factor procedure from scratch. In this special case, another approach called the method of undetermined coefficients is more effective. As we have seen in Theorem 3.5, the difference of two solutions of the nonhomogeneous equation (3.5) is a solution of the homogeneous equation (3.1). In turn, the sum of a solution of (3.5) and a solution of (3.1) is again a solution of (3.5). So, it is enough to find a single solution of (3.5) and the general solution of (3.1) and to add them in order to get the general solution of (3.5). In this context the single solution of (3.5) is often referred to as a particular solution, while the general solution of (3.1) is called the complementary solution. This so-called principle of superposition generally holds for linear equations (and only for linear equations). The method of undetermined coefficients looks for a special solution of (3.5) which has the same terms as I (t) and its derivatives. So, we try M(t) = A + B cos t + C sin t. UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES 23 A, B, C are the undetermined coefficients which give this method its name. Substituting this "ansatz" into (3.5) yields −B sin t + C cos t = −2(A + B cos t + C sin t) + 1 − cos t. Comparing the constant terms and the terms associated with sin t and cos t, we see that this equality holds if 0 = −2A + 1, −B = −2C, C = −2B − 1. So, A = 21 , C = − 15 , B = − 25 , and we have found a special solution of (3.5), M(t) = 1 2 − 2 5 cos t − 1 5 sin t. The general solution of (3.1) is given by N (t) = ke−2t and the general solution of (3.5) by N (t) = M(t) + ke−2t . We notice that N , while not approaching a constant, approaches the periodic function M as t → ∞. This is a phenomenon that holds in general for our scenario. Theorem 3.7. Consider N = ρ(t)N + I (t), N (t0 ) = N0 and assume that ρ is periodic with period p, ρ̄ = 0. (a) If ρ̄ < 0 and I is bounded on the interval [t0 , ∞), so is N . (b) If I is periodic, with same period p, then there exists some p-periodic function N∞ such that q(t) N (t) = [N (t0 ) − N∞ (t0 )]e(t−t0 )ρ̄ + N∞ (t), q(t0 ) p I (t − s) q(t) N∞ (t) = eρ̄s ds, ρ̄p 1−e q(t − s) 0 with q being the p-periodic function from the Floquet representation of the fundamental solution. For ρ̄ < 0, N (t) − N∞ (t) → 0 as t → ∞. If I is nonnegative and not identically 0, then N∞ is bounded away from 0 and has the opposite sign of ρ. So, if a fish pond in which the fish population cannot survive on its own is stocked periodically, the population size will eventually vary periodically as well (Figure 3.4C). In this context N∞ is called the persistent solution and N − N∞ the transient solution of the nonhomogeneous linear differential equation. However, if the population already grows on its own, then the stocking adds to the growth, which is exponential with periodic modulation (Figure 3.4B). 24 CHAPTER 3 Proof. (a) From the variation of constants formula (3.6) and Remark 3.2, t q(t) I (s) e(t−s)ρ̄ + q(t) ds. N (t) = N (t0 )e(t−t0 )ρ̄ q(t0 ) q(s) t0 Since q is periodic and strictly positive, it is bounded and bounded away from 0 on R. So, there exists a constant c > 0 such that q(t) c q(s) for all t, s ∈ R. So, for t t0 , N (t) cN (t0 ) + c t e(t−s)ρ̄ I (s) ds. t0 Furthermore, since I is bounded on [t0 , ∞), there exists some k > 0 such that I (s) k for all s t0 , thus 1 N (t) cN (t0 ) + ck (e(t−t0 )ρ̄ − 1). ρ̄ Since ρ̄ < 0, N (t) cN (t0 ) − ck . ρ̄ (b) We make the "ansatz" N∞ (t) = α t t−p Q(t) I (s) ds. Q(s) = ρ(t)N + I (t) and show We will determine α such that N∞ satisfies N∞ ∞ that N∞ has the form of part (b) of our theorem. Part (b) then follows from the fact that Ñ = N − N∞ satisfies Ñ = ρ(t)Ñ and from Theorem 3.5. Using the product rule (remember Q = ρ(t)Q) and the fundamental theorem of calculus, we see that Q(t) N∞ I (t − p) . (t) = ρ(t)N∞ (t) + α I (t) − Q(t − p) By Remark 3.2 and the periodicity of I , (t) N∞ eρ̄t q(t) = ρ(t)N∞ (t) + αI (t) 1 − ρ̄(t−p) . e q(t − p) Since q is p-periodic, (t) = ρ(t)N∞ (t) + αI (t)(1 − eρ̄p ). N∞ UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES 25 So, we choose α = (1 − eρ̄p )−1 . Notice that α has the opposite sign of ρ̄. Using Remark 3.2 in our ansatz, t q(t) N∞ (t) = α e(t−s)ρ̄ I (s) ds. q(s) t−p Substituting t − s = r we obtain the form in Theorem 3.7 (b), from which we realize that N∞ is p-periodic. Furthermore, we have the estimate N∞ (t) q(t)eρ̄p w(t) α with w(t) = t t−p I (s) ds. q(s) By the fundamental theorem of calculus, w (t) = I (t) I (t − p) − = 0, q(t) q(t − p) because I and q are p-periodic. This implies that w is constant and p I (s) w(t) = w(p) = ds > 0. 0 q(s) Since q is strictly positive (Remark 3.2), so is N∞ /α. More information can be obtained if ρ = ρ̄ = −α is constant (and negative). Then we can rewrite (3.5) as N = α(−N (t) + f (t)), where f is a periodic nonnegative function that is not identically 0, and use Fourier series. After scaling time we can assume that the period of f is 2π . If f is twice continuously differentiable, we can represent f by a Fourier series ∞ f (t) = ak eikt k=−∞ with ∞ |ak | < ∞. k=−∞ Recall that 1 a0 = 2π 0 2π f (t) dt = f¯ 26 CHAPTER 3 is the time average of f and the other Fourier coefficients are given by 2π 1 aj = f (t)e−ij t dt. 2π 0 We look for N = N1 + N2 with N1 being a general solution of N = −αN and N2 a particular solution in the form ∞ N2 = ak M k , k=−∞ Mk = α(−Mk + eikt ). We make the following "ansatz" for Mk , Mk = γk eik(t−tk ) . Substituting the ansatz into the differential equation yields after simplification that γk ike−iktk = −αγk e−iktk + α. After reorganization, α + ik γk . α The left-hand side has absolute value 1, so eiktk = γk = √ α α2 + k2 and α + ik eiktk = √ . α2 + k2 Separating real and imaginary parts (recall eiθ = cos θ + i sin θ ), tan(ktk ) = k , α which, for every k, has a unique solution tk ∈ [0, π/(2k)). Moreover, for every k = 0, tk → 0, α → ∞, tk → π/(2k), α → 0. Since N2 (t) = ∞ k=−∞ ak √ α α2 + k2 eik(t−tk ) and also the sum of the respective derivatives converge absolutely and uniformly, N2 is a solution of the differential equation. We see that the higher-frequency terms of N2 as compared with those of f are damped and that every term except UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES 27 the constant term for k = 0 lags behind that of f with a delay between 0 and a quarter of a period. The delays decrease with k. Let us consider the case where α is large. Since tk ≈ 0 except for large k, N2 ≈ f. Since the complementary solution is of the form N1 (t) = ce−αt , the population size N tracks the immigration rate. We turn to the case where α > 0 is small. As f is nonnegative and not identically 0, a0 > 0 and ak N2 (t) = a0 + α eik(t−tk ) . √ 2 2 α + k k=0 For small α, since ktk ≈ π/2 except for large k, ak N2 (t) ≈ a0 + α ei(kt−π/2) . k k=0 In this case, the population size exhibits small-amplitude oscillations of the order of α about its mean value with the fundamental mode of oscillation (k = 1) lagging behind that of the immigration rate by about one-quarter of a period. The population cannot react fast enough for its size to follow the oscillations of the immigration rate. Also t 1 ak ikt ¯ ¯ [f (s) − f ] ds − c N2 ≈ a0 + α e =f +α i k 0 k=0 for some constant c. This follows from the fact that 1 ak ikt e i k k=0 is an anti-derivative for f − f¯. To determine c, we observe from the differential equation for N2 , N2 = α(−N2 + f ), and the periodicity of N2 , that N2 and f have identical time averages. This implies that 2π t 1 [f (s) − f¯] ds dt. c= 2π 0 0 Changing the order of integration 2π 1 c= (2π − s)[f (s) − f¯] ds. 2π 0 Recalling the definition of f¯, 2π 1 ¯ c = πf − sf (s) ds. 2π 0 28 CHAPTER 3 Substituting this back into the equation for N2 , t 2π 1 N2 ≈ f¯ + α f (s) ds + (π − t)f¯ − sf (s) ds . 2π 0 0 3.2.2 Zero Average Intrinsic Growth Rate Let us briefly turn to the case where ρ is periodic with ρ̄ = 0 (i.e., on average over time, the per capita birth rate balances the per capita death and emigration rates) and that the immigration rate is periodic with the same period. Then t Q(t) I (s) + Q(t) ds. N (t) = N0 Q(t0 ) t0 Q(s) By Theorem 3.1 and Remark 3.2, Q is a periodic function with period p and so is the product I /Q. By Proposition 3.3, 1 p I (s) 1 t I (s) ds → ds =: κ, t → ∞. t t0 Q(s) p 0 Q(s) So, N (t)/(tQ(t)) → κ as t → ∞, in more suggestive formalism N (t) ≈ κtQ(t), t → ∞. This means that, for large times, the population size exhibits a linear increase in time which is modulated by a periodic function. Let us summarize the impact of immigration under periodicity. For an illustration, compare Figures 3.3 and 3.4, where we use the normalized per capita birth and death rates, and the immigration rate in Figure 3.2. If the average intrinsic growth rate ρ̄ is negative (case C in Figures 3.3 and 3.4), immigration turns a population that would become extinct otherwise into an eventually periodic population. If the time average ρ̄ is zero (case A), immigration converts a periodic population into one that grows linearly with periodic modulation. If ρ̄ > 0 (case B), the population size increases exponentially with exponent ρ̄ and periodic modulation, independently of whether or not there is immigration, but immigration adds to the growth. Exercises 3.1. Consider unrestricted population growth in a closed population with constant per capita birth and death rates. Calculate the doubling time of the population in terms of ρ̄ if ρ̄ > 0. Remark. The doubling time T is characterized by N (t0 + T ) = 2N (t0 ). 29 UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES A B 180 160 140 120 100 80 60 40 20 0 30 25 y(t) 20 15 10 5 0 0 2 4 t 6 8 10 0 2 4 6 8 10 t 6 8 10 C 4 y(t) 3 2 1 0 2 4 t Figure 3.4. Dynamics of open populations under various combinations of the per capita birth and death rates of Figure 3.2 and the immigration rate of Figure 3.2. The averages of the intrinsic growth constant have been chosen as in Figure 3.3. (A) ρ̄ = 0. (B) ρ̄ = 0.3. (C) ρ̄ = −0.9. 3.2. If the per capita mortality rate µ is constant, then L = 1/µ is the life expectancy (expected age at death) of the individuals in the population (see (2.3)). If the per capita birth rate β is also constant, then R0 := βL = β/µ is the average number of offspring that one typical individual produces during its lifetime. R0 is called the basic reproduction ratio of the population. In a population with two sexes, one would divide the average number of offspring by two or only count the average number of daughters per female (in other words one would only consider the female part of the population). Express the doubling time of the population in terms of R0 and L. Calculate the doubling time if the life expectancy L is 70 years and R0 = 1.5 (i.e., in a human population, every female would have an average number of 1.5 daughters during her lifetime). 3.3. Determine the solution of N = tan(t)N, N (0) = N0 . 30 CHAPTER 3 Where is the solution defined? How does it behave? Where does it make biological sense? 3.4. (a) Solve N = −αN + A with positive constants α and A and initial data N (t0 ) = N0 . (b) Using (a), show that N (t) → A/α as t → ∞. 3.5. Determine the solution of N = ρN + I (t), N (0) = N0 with constant ρ and immigration rate I (t) = 1 + sin t. What can you say about the behavior for large times? 3.6. Consider N = sin(t)N + 1, N (0) = N0 . What can you tell about the growth of the population as t → ∞ from the theory developed in the text. Try the problem on a computer for various initial values N0 . 3.7. (a) Solve N = a N + bt γ , t > 0, t N (t0 ) = N0 , with constants a ∈ R, γ ∈ R, b > 0 and initial data N0 0, t0 > 0. (b) Discuss the behavior of the solution as t → ∞. 3.8. Solve N = ρ(t)N + I (t), N (t0 ) = N0 for ρ(t) = − sin t − 1, I (t) = a(1 + sin t) with some constant a > 0. Discuss the behavior as t → ∞. 3.9. Prove the second statement in Theorem 3.4 that the size of a closed population decreases exponentially if the asymptotic time average of ρ exists and ρ̄ < 0. 3.10. A continuous function f : R → R is called quasi-periodic if it is the finite sum of periodic functions with different periods. Show that the asymptotic time average exists for quasi-periodic functions. UNCONSTRAINED POPULATION GROWTH FOR SINGLE SPECIES 31 3.11. A continuous function f : [0, ∞) → [0, ∞) is called asymptotically almost periodic if for every > 0 there exists a quasi-periodic function g and some t > 0 such that |f (s) − g(s)| , for all s t. Show that the asymptotic time average exists for asymptotically almostperiodic functions. Hint: consider the sequence 1 xn = n n f (s) ds 0 and show that it is a Cauchy sequence. Chapter Four Von Bertalanffy Growth of Body Size Here we do not consider the growth of a population, but the size growth of a single individual. We include this topic for two reasons. On the one hand, we can consider the body to be a population of cells so that this topic fits under the umbrella of population growth. On the other hand, growth models for the body sizes of individuals are the basis for size-structured population models, which stratify a population along the body sizes of its constituents (Metz, Diekmann, 1986; Kooijman, 1993, 2000). The von Bertalanffy model (1934) is built on the assumption that the acquisition of food is proportional to the surface of the body, S, while maintenance costs for the body are proportional to the volume of the body, V : V = F (t)S − ν(t)V . F (t) is the quantity of food taken up per unit surface area and ν(t) denotes the maintenance costs per unit volume. Since S = γ V 2/3 with a proportionality constant γ that depends on the shape of the body, we have (4.1) V = γ F (t)V 2/3 − ν(t)V . Since γ does not depend on time, we have implicitly assumed that, though growing, the body does not change its shape. In order to transform this equation into a linear one, we introduce the length of the body, L. Since V = (αL)3 with some proportionality constant α depending on the shape of the body, we have V = 3α 3 L2 L . Thus we obtain the following ordinary differential equation for L: 3α 3 L2 L = γ F (t)α 2 L2 − ν(t)α 3 L3 . Hence L = f (t) − µ(t)L, (4.2) with f (t) = (γ /3)α −1 F (t), µ(t) = ν(t)/3. Equation (4.2) can now be solved by formula (3.6). In a laboratory environment, one can assume that the organism under question is fed with a period of one day. Furthermore, since the maintenance costs often depend on temperature, µ may have the same period. From Theorem 3.7 and V = (αL)3 we immediately obtain the following result. 34 CHAPTER 4 Theorem 4.1. Let F and ν be p-periodic nonnegative continuous functions that are not identically 0. Then there exists a strictly positive periodic function V∞ such that V (t) − V∞ (t) → 0 as t → ∞. Let F and ν, and so f and µ, be constant. If L0 is the length at time t = 0, the length at time t is given by L(t) = L0 e−µt + L∞ (1 − e−µt ) = (L0 − L∞ )e−µt + L∞ , (4.3) f L∞ = . µ This shows that the length converges to L∞ = f/µ as t → ∞, which can be interpreted as the final length under constant conditions. If we go backwards in time and if L0 < L∞ , we notice that L becomes negative at t− = 1 L∞ − L0 , ln µ L∞ and the solution stops making sense. If we go back to body volume, however, we find that the following is also a solution to (4.1) which is defined for all times: t t− , (L0 e−µt + L∞ (1 − e−µt ))3 , (4.4) V (t) = 0, t t− . This is, of course, related to the fact that (4.1) can have more than one solution in a neighborhood of V = 0, because the vector field fails to be Lipschitz there. While the graph of the solution L in (4.3) is concave, the graph of V in (4.4) is sigmoid. We will see more of this later. The von Bertalanffy model is an example of modeling from first principles. The principles are that the per capita food intake scales with body surface, while the maintenance costs scale with body volume. For small creatures it may be difficult to verify these principles directly, and it is even possible that they may not tell the whole story. So, it is imperative to check whether the above formulas provide results that agree with real data of body sizes. To this end one needs to know the parameters in (4.3). We now interpret t as age and L0 as length at birth, which occurs at t = 0. By observing for long enough, it is possible to get a good idea of the final body length L∞ , while it may be difficult to determine the length at birth (the individuals may yet be too fragile to be measured). The parameter µ, which is related to maintenance costs, is very difficult to determine directly. The idea is to find L0 and µ such that the length growth data and the outputs of equation (4.3) deviate from each other as little as possible. Rewriting (4.3) as L(t) − L∞ = (L0 − L∞ )e−µt , 35 VON BERTALANFFY GROWTH OF BODY SIZE we can take logarithms, ln(L∞ − L(t)) = ln(L∞ − L0 ) − µt. The right-hand side is a line with slope −µ. Set b(t) = ln(L∞ − L(t)), b0 = ln(L∞ − L0 ) and assume that we have n measurements b(tj ) = bj at ages t1 , . . . , tn . This gives us the linear system bj = b0 − µtj for j = 1, . . . , n, which takes the matrix form Ax = b with b = (b1 , . . . , bn )T , x = (b0 , −µ)T and 1 ··· 1 . AT = t1 · · · tn The superscript "T" denotes the transpose of a matrix. For n > 2 this is an overdetermined system which typically has no solution. However, it can be shown by projection methods (Edwards, Penny, 1988, Sections 5.2, 5.3, for example) that the system AT Ax ∗ = AT b has a solution x ∗ and that the Euclidean distance Ax ∗ − b is the minimum of all distances Ax − b , x ∈ R2 . x ∗ is called the least-squares solution of Ax = b and the respective procedure least-squares fitting. Once we have extracted L0 and µ from x ∗ , we can solve L = µ(L∞ − L), L(0) = L0 numerically and see whether the plotted solution graph comes sufficiently close to the measurements. Typically, the von Bertalanffy model fits individual growth curves under constant food regimes quite well. See Metz, Diekmann (1986, A.I.3.1), Kooijman, Metz (1984) and Kooijman (1986). For refinements of this model and more growth curves see the book by Kooijman (1993, 2000). Exercises 4.1. Determine whether the following growth data of body length resemble von Bertalanffy growth under constant food regime: age (days) length (mm) 4 8 9.5 12 16 1.8 2.3 2.5 2.6 2.8 The final body length is assumed to be L∞ = 3 mm. 4.2 Consider V = ν(t)(AV 2/3 − V ) with a positive constant A and a nonnegative, continuous function ν defined on R. Assume that the asymptotic time average ν̄ exists and ν̄ > 0. What will be the final body volume (with proof )? Chapter Five Classic Models of Density-Dependent Population Growth for Single Species For many natural populations, the available resources are limited and, directly or indirectly, the size of the population couples back to the per capita birth and/or mortality rates. The simplest way to incorporate this feedback consists in assuming that the per capita birth and mortality rates at time t depend on the population size at that very moment t: B(t) = β(t, N (t))N (t), D(t) = µ(t, N (t))N (t). For a closed population, equation (2.1), with migration E = I = 0, takes the form N = N [β(t, N ) − µ(t, N )]. (5.1) We will present the classical equations by Bernoulli and Verhulst (the second being a special case of the first), by Ricker, Beverton–Holt, Smith, and Gompertz, and analyze them as special cases of a certain class of nonlinear scalar ordinary differential equations. We will treat them in some detail because they not only model the one-species situation, but often occur also as modules in multi-species models. 5.1 The Bernoulli and the Verhulst Equations A first class of models assumes that the population density only affects the per capita mortality rate, µ(t, N) = m0 (t) + m(t)N θ (5.2) where θ is a positive constant. Here m0 is the density-independent component of the per capita mortality rate. Then N = ρ(t)N − m(t)N θ+1 (5.3) ρ(t) = β(t) − m0 (t). (5.4) with This is (Jakob) Bernoulli's equation. Leibnitz (1696) suggested the following transformation to reduce it to a linear equation (see the problems in Section 2.2 38 CHAPTER 5 of Boyce, DiPrima, 1992): x = N −θ . Then x = −θ N −θ −1 N = −θN −θ−1 [ρ(t)N − m(t)N θ+1 ] = −θρ(t)x + θ m(t). This equation is of the form (3.5), and a solution x can be found from formula (3.6). In turn, N = x −1/θ is a solution of (5.4). Again we may like to model seasonal influences on the birth and mortality rates by periodic functions (with the period being one year, for example). Then we obtain the following result from Theorem 3.7 (b). Theorem 5.1. Let θ > 0 and let ρ(t), m(t) be continuous p-periodic functions. Assume that 1 p ρ(t) dt > 0 ρ̄ = p 0 and that m is nonnegative and not identically 0. Then there exists a strictly positive p-periodic function N∞ such that N (t) − N∞ (t) → 0 as t → ∞ for every solution N of (5.3). If ρ and m are time independent and positive, it is convenient to rewrite (5.3) as N = ρN (1 − (N/K)θ ) (5.5) K = (ρ/m)1/θ . (5.6) with K is called the carrying capacity for the population, because N (t) → K as t → ∞. (This follows from Corollary 3.6 applied to x.) In Chapter 6 we will see that this is also related to the fact that K is the unique positive zero of the right-hand side of (5.5). For the case θ = 1, equation (5.5) has become known as Verhulst's equation, the Verhulst–Pearl equation, or the logistic equation as Verhulst (1845) called it himself. Gilpin and Ayala (1973) have rediscovered (5.5) for θ = 1, and theoretical ecologists sometimes refer to it as the asymmetric or generalized logistic equation rather than Bernoulli's equation (Hutchinson, 1978, Chapter 1). As we will see later in a more general context, the solution graph of the Verhulst equation, for any initial value between 0 and K, has an S-shape or sigmoid shape, and connects 0 (for t → −∞) with K (for t → ∞), and has an inflection point at 21 K. But the symmetry of the solution graph with respect to 1 2 K (half the carrying capacity) is even more striking. Proposition 5.2. Let N be the solution of the Verhulst equation N = ρN (1 − N/K) with ρ, K > 0. Without restriction of generality let N (0) = 21 K. Then N (t) − 21 K = 21 K − N (−t). CLASSIC MODELS OF DENSITY-DEPENDENT POPULATION GROWTH 39 Proof. We set x(t) = N (t) − 21 K, y(t) = 21 K − N (−t). Obviously, x(0) = 0 = y(0). We shall show that x and y satisfy the same differential equation with a locally Lipschitz continuous vector field. Then x = y by the uniqueness theorem for ordinary differential equations (see Hale (1980), for example). Indeed, x + 21 K K 1 x 1 x = N = ρ(x + 2 K) 1 − =ρ x+ − , K 2 2 K 1 K − y(t) K y (t) = N (−t) = ρ − y(t) 1 − 2 2 K 1 y K =ρ − +y . 2 K 2 It was this striking symmetry that led Gilpin and Ayala (1973) to consider the asymmetric logistic alias Bernoulli equation. Let us explore the symmetry properties of solutions N to (5.5) for θ > 0. Set z = N θ . Then z z = θ N θ−1 N = θρz 1 − θ . K Hence z is a solution of the Verhulst equation with carrying capacity K θ . By Proposition 5.2, if z(0) = K θ /2, z(t) − K θ /2 = K θ /2 − z(−t). Going back to N we have the following result. Proposition 5.3. Let N be a solution to (5.5) with ρ, K, θ > 0. Let N (0) = K2−1/θ . Then (N (t))θ + (N (−t))θ = K θ ∀t 0. 5.2 The Beverton–Holt and Smith Differential Equation Other models assume that the population density affects the (effective) per capita birth rate rather than the per capita mortality rate. The effective per capita birth rate takes account of higher than usual mortality at or immediately after birth. For fishery models, though in a somewhat different context, Ricker (1954, 1975) suggested β(N ) = β0 e−αN , 40 CHAPTER 5 while Beverton, Holt (1957) suggested β(N ) = β0 . 1 + αN As we will see, both functional forms can be explained by heavy cannibalization of juveniles by adults of the same species, during a very short period immediately after birth. The different forms originate from different assumptions concerning the length distribution of the juvenile period; in the Ricker equation the juvenile period has a fixed length, while in the Beverton–Holt equation its length is exponentially distributed, i.e., individuals make the transition from juvenile to adult stage at a fixed constant rate. In a way, it is disconcerting that different assumptions about the length distribution of the juvenile period lead to these dramatically different forms (see Section 13.7 for an elaboration of this theme). Other explanations than cannibalism are possible, in particular for the Beverton–Holt reproduction function, which can also be obtained from a resource–consumer model via a time scale argument. We will deal with this explanation first and turn to cannibalism later. 5.2.1 Derivation from a Resource–Consumer Model Let M be the biomass of the population under consideration (the consumers), and F the biomass of the food (resource) on which the population lives: F = Λ − νF − b F M, γ M = bF M − µM. The food is provided at the constant rate Λ and degrades at a constant per unit rate ν. Food consumption follows the law of mass action, i.e., the consumption rate is proportional to both the biomass of food available and the biomass of consumers; the proportionality constant b gives the per unit consumer biomass gained from one unit biomass of food per unit of time. The yield constant γ describes the conversion of consumed food biomass into consumer biomass, i.e., γ is the per unit consumer biomass gained from one unit of consumed food. The quotient b/γ is the per unit food biomass consumed per unit of consumer biomass and per unit of time. We mention that the biomass gain of the consumers is not only due to the increase of individuals in weight or volume, but also to the production of offspring. Assuming that the food dynamics are much faster than the consumer dynamics, we want to make a quasi-steady-state approach to the food equation. In order to be able to compare the speeds of the different dynamics we introduce new, dimensionless variables. We start with the dependent variables. Without CLASSIC MODELS OF DENSITY-DEPENDENT POPULATION GROWTH 41 the consuming population being around, the food would converge to Λ/ν, so we introduce Λ F = z. ν Substituting this expression into the differential equations, Λ bΛ z = Λ − Λz − zM, ν γ ν bΛ M = µM z−1 . µν Simplifying and introducing the basic biomass production ratio, R0 = yields Λb , ν µ (5.7) b z =ν 1−z− zM , γν M = µM(R0 z − 1). The form of R0 motivates us to define M= Λγ x, µ and obtain z = ν(1 − z − R0 zx), x = µx(R0 z − 1). To interpret R0 , recall that Λ/ν is the large-time limit of the food concentration if there are no consumers around. Since b is the per unit/unit biomass gain rate and 1/µ is the life expectancy of the consumer, we see that R0 is the ratio by which consumer biomass grows over a consumer's lifetime, without competition from conspecifics. Finally, we introduce the dimensionless time τ , τ = µt; in this way our time unit becomes the life expectancy of the consumer, 1/µ. Let ż and ẋ denote the derivatives with respect to the new time τ , by the chain rule ż = µz and ẋ = µz . Dividing the first equation by ν and the second by µ we obtain ż = 1 − z − R0 zx, ẋ = x(R0 z − 1), with = µ/ν. 42 CHAPTER 5 The assumption that the food dynamics are fast compared to the population dynamics can now be rephrased more precisely as 1 and R0 being not too small compared to 1. We expect this system to behave similarly to the system 0 = 1 − z − R0 zx, ẋ = x(R0 z − 1). This formal argument can be made rigorous using singular perturbation theory (Tikhonov et al., 1985; Hoppensteadt, 1974). From the first equation we get the quasi-steady-state 1 z= , 1 + R0 x which we substitute into the second equation: R0 ẋ = x −1 . 1 + R0 x (5.8) We return to unscaled time and consumer biomass and obtain the Beverton– Holt equation, R0 M = µM −1 1 + (R0 µ/Λγ )M β =M −µ , (5.9) 1 + αM with β = µR0 = Λ b, ν α= β . Λγ 5.2.2 Derivation from Cannibalism of Juveniles by Adults In order to derive the Beverton–Holt reproduction function from adults cannibalizing juveniles, let N (t) denote the size of the adult population and J (t) the size of the juvenile population. Then J = β0 N − ν0 J − κJ N − γ J, N = γ J − µN. Here β0 is the per capita birth rate of the population, ν0 and µ are the natural per capita mortality rates of juveniles and adults, respectively, γ is the per capita rate of a juvenile turning into an adult, and κ is the per capita rate of juveniles being cannibalized by adults per capita of adults per unit of time. The cannibalism rate is modeled according to the law of mass action, i.e., it is proportional to both the size of the adult and the juvenile populations. 43 CLASSIC MODELS OF DENSITY-DEPENDENT POPULATION GROWTH 1/γ is the average length of the juvenile period (death disregarded), and we assume that it is very short compared to the life expectancy of an adult, 1/µ, i.e., µ = γ , with > 0 being a very small dimensionless parameter. We also assume that cannibalism occurs on a similar time scale to leaving the juvenile stage, i.e., κ = α̃ for some parameter α̃. Let J = xN . Then x + xµ(x − 1) = β0 − ν0 x − α̃xN − µx, N = µN (x − 1). Finally, we introduce the dimensionless time τ , τ = µt; in this way our time unit becomes the life expectancy of the adult population, 1/µ. Let ż and ẋ denote the derivatives with respect to the new time τ , by the chain rule ż = µz and ẋ = µz . Dividing both equations by µ we obtain ẋ + x(x − 1) = β0 ν0 α̃ − x − xN − x, µ µ µ Ṅ = N (x − 1). Since is very small, we assume the solution of the first equation to behave as if = 0, 0 = β0 − α̃xN − µx. Solving this equation, we obtain the quasi-steady-state x= β0 , µ + α̃N the substitution of which into the equation N = µ(x − 1) again yields equation (5.9) with α = α̃/µ. If β > µ, we can give equation (5.9) the following form suggested by Smith (1963) (via a different derivation (see also Hutchinson, 1978, Chapter 1)): M = ρM 1 − (M/K) 1 + αM with ρ = β − µ, K= (5.10) β −µ . αµ Notice that we formally recover the Verhulst equation (also known as the logistic equation) from (5.10) setting α = 0. We realize that all positive solutions remain bounded, since M < 0 for M > K; in fact M(t) max{M(0), K} for all t 0. 44 CHAPTER 5 Let us try to find a solution of the Beverton–Holt equation in the form of Smith's equation. Equation (5.10) has two constant solutions, M ≡ 0 and M ≡ K. To obtain other solutions we separate the variables, 1 + αM 1 − (M/K) + (α + (1/K))M = M M[1 − (M/K)] M[1 − (M/K)] M 1 M (M/K) (M/K) = + α+ = + (1 + αK) . M K 1 − (M/K) (M/K) 1 − (M/K) ρ = M Integrating, |M| M − (1 + αK) ln1 − , (5.11) K K with an integration constant k. It seems impossible to find an explicit solution unless α = 0, in which case Smith's equation reduces to the logistic equation. But we see from this equation that M(t) → K as t → ∞, because this is the only way in which the right-hand side can approach infinity (recall that all positive solutions are bounded). So, as with the logistic equation, K plays the role of the carrying capacity. We can also use relation (5.11) to test whether given population data can be fitted by Smith's equation, provided that we can observe the population for long enough to get a good idea what the carrying capacity K may be. Let x = (k, ρ, 1 + αK)T be the vector that contains the parameters we want to determine and Mj = M(tj ), j = 1, . . . , n, be measurements taken at times t1 , . . . , tn . Then we have the linear system ρt + k = ln x1 + tj x2 + ln |1 − (Mj /K)|x3 = ln |Mj /K|, j = 1, . . . , n, in matrix form Ax = b with bT = (ln |M1 /K|, . . . , ln |Mn /K|), 1 ··· 1 . AT = t1 ··· tn ln |1 − (M1 /K)| · · · ln |1 − (Mn /K)| If n > 3, this is an overdetermined system which typically has no solution; its least-squares solution can be found by solving the linear system AT Ax = AT b (cf. Chapter 4). Fortunately this has already been programmed into packages like Maple, Mathematica or Excel. If we restrict our fitting exercise to the Verhulst equation, α = 0, we obtain the relation Mj , j = 1, . . . , n. x1 + tj x2 = ln K − Mj By (5.11), K − M always has the same sign; most data are collected for the case in which M takes values between 0 and K, x1 + tj x2 = ln Mj , K − Mj j = 1, . . . , n. CLASSIC MODELS OF DENSITY-DEPENDENT POPULATION GROWTH 45 If the data are collected at equidistant times tj = t1 +(j −1)h, we can try to use the following procedure for estimating K, which has already been suggested by Verhulst (1845, Section 8). By subtracting the last equation from the analogous one for j + 1, we obtain hx2 = ln Mj +1 Mj − ln , K − Mj +1 K − Mj j = 1, . . . , n − 1. Equating the right-hand sides and rearranging, ln Mj +2 Mj Mj +1 + ln = 2 ln , K − Mj +2 K − Mj K − Mj +1 j = 1, . . . , n − 2. We exponentiate and rearrange, Mj2+1 (K − Mj +2 )(K − Mj ) = Mj +2 Mj (K − Mj +1 )2 . We expand this expression and solve for K, K = Mj +1 Mj +2 Mj +1 + Mj +1 Mj − 2Mj +2 Mj , Mj2+1 − Mj +2 Mj j = 1, . . . , n − 2. For real data, we will obtain different values for K as j = 1, . . . , n − 2, so we could take their average. If a population grows less than exponentially, one has Mj +2 Mj +1 < , Mj +1 Mj i.e., the denominator in the expression for K is positive. If the data are taken from a regime, however, where the growth is still close to being exponential, then the denominator is close to 0. The numerator can be rewritten as Mj +2 Mj Mj +2 + Mj Mj +1 − Mj +2 Mj . 2Mj +1 − Mj +2 Mj + Mj +1 2 Since the arithmetic mean is larger than the geometric mean, this expression is positive and safely bounded away from 0. This means that our formula for K is very sensitive to errors in measurement and hardly practical, as long as the population growth is still almost exponential. In fact, it can easily happen that one obtains negative values if one uses real data. But even for data which level off for large times, I have found this procedure impractical because the values for K that one obtains vary too widely. 5.3 The Ricker Differential Equation While we have assumed before that juveniles turn into adults at a constant rate, we now assume that the juvenile stage has a fixed length, τ , which is very short 46 CHAPTER 5 compared to the life expectancy of adults, 1/µ. So, individuals that enter the adult population at time t have been born at time t − τ , and we obtain the equation N (t) = β0 N (t − τ )P (t, t − τ ) − µN. As before, β0 is the per capita birth rate and µ the per capita mortality rate of adults. P (t, r) denotes the probability of a juvenile to survive from time r to time t, t r. A similar derivation as for (2.2) expresses P in terms of a per capita juvenile mortality rate ν, t P (t, r) = exp − ν(s) ds . r We assume that the per capita juvenile mortality rate splits into a natural mortality rate and a mortality rate due to cannibalism, ν(t) = ν0 + κN (t). Here ν0 is the natural per capita mortality rate, and κ is the per capita rate at which juveniles are cannibalized by adults, per capita of adults. In other words, 1/κ is the average time it takes a juvenile to be killed by one unit of adults. As in the previous section, we assume that 1/κ is on the same time scale as the length of the juvenile period, κ = α/τ. By substitution into the formula for P , P (t, t − τ ) = e −ν0 τ α exp − τ t N (s) ds . t−τ Taking the limit for τ → 0, P (t, t−) = e−αN(t) . Taking also the limit τ → 0 in the original differential equation for N , we obtain the Ricker equation N = β0 N e−αN − µN. (5.12) Let us rewrite (5.12) in a form that displays the carrying capacity. Factoring out µ, β0 −αN N = µN −1 . e µ Set 1 β0 K = ln . α µ Then N = µN (eαK(1−(N/K)) − 1). CLASSIC MODELS OF DENSITY-DEPENDENT POPULATION GROWTH 47 Now set ρ = µ(eαK − 1), Then γ = αK. ρ = β0 − µ, and N = ρN β0 γ = ln , µ eγ (1−(N/K)) − 1 . eγ − 1 Notice that we recover the logistic equation as the limit of this equation for γ → 0 (apply l'Hôpital's rule). K again plays the role of the carrying capacity as the right-hand side of this equation is 0 for N = K. 5.4 The Gompertz Equation The equation N = rN ln(K/N ), with positive constants r, K > 0, which ultimately goes back to Gompertz (1825), does not fit into the fundamental balance equation of population dynamics, equation (2.1), with per capita birth and mortality rates. It has been very successful in fitting data of tumor growth, with N representing the number of tumor cells. This equation can be solved explicitly. Let x = ln(N/K) = ln N − ln K. Then x = N /N = −rx. Thus ln |x| = −rt + k and ln(|ln(N/K)|) = k − rt. From this equation we see that N (t) → K as t → ∞ and N (t) → 0 as t → −∞ if N (0) < K, while N (t) → ∞ as t → −∞ if N (0) > K. As in the logistic equation, we can estimate K from the large-time behavior of the data and r and k by least-squares fitting as in Chapter 4. This is a further reason for the popularity of the Gompertz equation. 5.5 A First Comparison of the Various Equations In order to compare these various differential equations in the time-autonomous case, we write them in the respective forms that display the carrying capacity, 48 CHAPTER 5 A 1.0 B 0.8 x V V B+ B− Sa Sb 0.6 0.4 0.2 0 C 1.0 D V 0.8 x 0.6 G Ra Rb V 0.4 0.2 0 −4 −2 0 time 2 4 −4 −2 0 time 2 4 Figure 5.1. Comparison of the Verhulst (logistic) equation (labelled "V" and shown in all parts of the figure), (A) the Bernoulli equation, (B) the Smith equation, (C) the Ricker equation, and (D) the Gompertz equation. (A) V, θ = 1; B−, θ = 0.5; B+, θ = 2. (B) V, α = 0; Sa, α = 1; Sb, α = 3. (C) V, γ → 0; Ra, γ = 1; Rb, γ = 3. See the text for further explanation. K. Scaling the population size, we can achieve that K = 1. Scaling time, we can make an additional parameter equal to 1, ρ = 1 or r = 1, respectively: x = x(1 − x) (Verhulst, logistic), x = x(1 − x θ ), θ >0 1 − x x = x , α0 1 + αx eγ (1−x) − 1 x = x , γ 0 eγ − 1 x = −x ln x (Bernoulli), (Smith, Beverton–Holt), (Ricker), (Gompertz). As we have noticed before, theVerhulst equation is a special case of the Bernoulli equation with θ = 1, and a limit case of the Smith equation with α = 0 and of the Ricker equation with γ → 0. In Figure 5.1 we have plotted solutions for the Bernoulli equation with θ = 1 (Verhulst), θ = 0.5, and θ = 2; for the Smith equation with α = 0 (Verhulst), α = 1, and α = 3; and for the Ricker equation with γ → 0 (Verhulst), γ = 1, and γ = 3. CLASSIC MODELS OF DENSITY-DEPENDENT POPULATION GROWTH 49 All solutions show an S-shaped growth, often called sigmoid growth. They increase from one equilibrium state, x = 0, to another one, x = 1, first in a convex and then in a concave manner. Notice that in all other examples the symmetry of the Verhulst equation is lost; in particular, the inflection point is no longer at x = 0.5. Only the Bernoulli equation is flexible enough to have the inflection point at both x > 0.5 (Figure 5.1A, B− for θ = 0.5) and x < 0.5 (Figure 5.1A, B+ for θ = 2). In the Smith, Ricker and Gompertz equations, the inflection point is at x < 0.5. In Chapter 6, we will find that these are general features that do not depend on the particular choices of the parameters we have made. See the remark between Lemma 6.2 and Theorem 6.4. Chapter Six Sigmoid Growth Rather than analyzing the many differential equation models we derived in Chapter 5 individually, we try a unified approach. First we notice that the qualitative behavior of solutions of Mathematics In Population Biology Thieme Pdf
Source: https://b-ok.asia/book/18153138/d24666
Posted by: griffinfloont.blogspot.com

0 Response to "Mathematics In Population Biology Thieme Pdf"
Post a Comment