Mathematical Modeling

Comprehensive guide to mathematical modeling covering model construction, analysis, validation, and applications across science, engineering, and economics.

mathematical-modeling differential-equations optimization probability-models systems-analysis applied-mathematics

Prerequisites

Required:

Calculus - Derivatives, integrals, and their applications
Differential Equations - Solution methods and applications
Algebra - Functions, systems of equations

Helpful:

Statistics - For stochastic models and data analysis
Numerical Methods - For computational solutions
Matrices - For linear models and systems
Domain knowledge in area of application

Overview

Mathematical modeling is the process of using mathematical language and tools to describe, analyze, and predict real-world phenomena. It bridges the gap between theoretical mathematics and practical applications, providing insights into complex systems across science, engineering, economics, biology, and social sciences.

This comprehensive guide covers the modeling process, various types of models, techniques for model construction and validation, and applications across multiple domains. Learn to translate real problems into mathematical frameworks and extract meaningful insights.

The Modeling Process

Stages of Mathematical Modeling

1. Problem Identification
   • Define the real-world problem clearly
   • Identify key variables and relationships
   • Determine modeling objectives

2. Model Formulation  
   • Make simplifying assumptions
   • Choose appropriate mathematical tools
   • Construct mathematical equations

3. Mathematical Analysis
   • Solve equations analytically if possible
   • Use numerical methods for complex models
   • Analyze model behavior and properties

4. Model Validation
   • Compare predictions with observed data
   • Test model sensitivity to parameters
   • Assess model limitations and assumptions

5. Implementation and Interpretation
   • Apply model to solve original problem
   • Interpret results in real-world context
   • Communicate findings to stakeholders

6. Model Refinement
   • Improve model based on validation results
   • Add complexity if needed
   • Iterate through process as necessary

Modeling Principles

Occam's Razor: Prefer simpler models that explain the data adequately
"All models are wrong, but some are useful" - George Box

Key Considerations:
• Purpose: What questions should the model answer?
• Scope: What phenomena should be included/excluded?
• Scale: Temporal and spatial scales of interest
• Accuracy: Required precision for intended use
• Tractability: Mathematical solvability
• Interpretability: Can results be understood and explained?

Trade-offs:
• Realism vs. Simplicity
• Accuracy vs. Tractability  
• Generality vs. Specificity
• Detail vs. Insight

Types of Mathematical Models

Deterministic vs. Stochastic Models

Deterministic Models:
• Same inputs always produce same outputs
• Use differential equations, algebraic equations
• Examples: Newton's laws, population growth models

Stochastic Models:
• Include randomness and uncertainty
• Use probability distributions, random variables
• Examples: Stock prices, disease spread with uncertainty

Example Comparison - Population Growth:
Deterministic: dP/dt = rP (exponential growth)
Stochastic: dP = rP dt + σP dW (adds random fluctuations)
where dW represents random "noise"

Discrete vs. Continuous Models

Discrete Models:
• Time/space divided into discrete steps
• Use difference equations, matrices
• Suitable for: populations counted in integers, time steps

Continuous Models:  
• Variables change smoothly over continuous domains
• Use differential equations, calculus
• Suitable for: physical quantities, large populations

Example - Compound Interest:
Discrete: A(n) = P(1 + r)^n (annual compounding)
Continuous: A(t) = Pe^(rt) (continuous compounding)

Linear vs. Nonlinear Models

Linear Models:
• Variables appear to first power only
• Superposition principle applies
• Often easier to solve analytically

Examples:
• Spring force: F = -kx
• Linear population growth: dP/dt = r
• Linear regression: y = mx + b

Nonlinear Models:
• Variables raised to powers, multiplied together, in transcendental functions
• Rich, complex behavior possible
• Often require numerical solutions

Examples:
• Logistic growth: dP/dt = rP(1 - P/K)
• Pendulum: d²θ/dt² = -(g/L)sin(θ)
• Predator-prey: dx/dt = ax - bxy, dy/dt = cxy - dy

Static vs. Dynamic Models

Static Models:
• Describe equilibrium or steady-state situations
• No time dependence
• Use algebraic equations

Examples:
• Supply and demand equilibrium
• Structural engineering stress analysis
• Optimization problems

Dynamic Models:
• Describe how systems change over time
• Include time as independent variable
• Use differential equations, difference equations

Examples:
• Population dynamics
• Chemical reaction kinetics  
• Economic growth models

Basic Modeling Techniques

Dimensional Analysis

Uses principle that physical equations must be dimensionally consistent

Process:
1. Identify relevant variables and their dimensions
2. Form dimensionless groups (π-terms)
3. Express relationships between dimensionless groups

Example: Pendulum Period
Variables: T (period), L (length), g (gravity), m (mass)
Dimensions: [T] = T, [L] = L, [g] = LT⁻², [m] = M

Dimensional analysis shows T = f(L/g) where f is dimensionless
Physical reasoning suggests T ∝ √(L/g)
Exact solution: T = 2π√(L/g)

Benefits:
• Reduces number of variables
• Guides experimental design
• Provides insight into scaling laws

Conservation Principles

Many models based on conservation laws:
• Conservation of mass (continuity equation)
• Conservation of energy  
• Conservation of momentum
• Conservation of charge
• Conservation of species

General Form: Rate of Change = Input - Output + Generation

Example: Population Model
Rate of change = Births - Deaths + Immigration - Emigration
dP/dt = bP - dP + I - E = (b - d)P + (I - E)

Example: Heat Transfer
dT/dt = (Heat in - Heat out + Heat generated) / (mass × specific heat)

Example: Chemical Reactor
dC/dt = (Flow in × C_in - Flow out × C) / Volume + Reaction rate

Scaling and Limiting Behavior

Understanding model behavior in extreme cases:

Small Parameter Analysis:
If ε << 1, expand solutions in powers of ε
Leading order behavior often provides key insights

Large Parameter Analysis:  
As parameter → ∞, what happens to solution?

Example: Damped Oscillator
m(d²x/dt²) + c(dx/dt) + kx = 0

Light damping (c small): oscillatory behavior
Heavy damping (c large): exponential decay
Critical damping: fastest return to equilibrium

Scaling reveals dimensionless groups:
ω₀ = √(k/m), ζ = c/(2√(km))
Behavior depends only on ζ, not individual parameters

Population Dynamics Models

Exponential Growth Model

Assumptions:
• Unlimited resources
• Constant birth and death rates
• No age structure

Model: dP/dt = rP
Solution: P(t) = P₀e^(rt)

Where:
• P(t) = population at time t
• P₀ = initial population  
• r = intrinsic growth rate (births - deaths per capita)

Properties:
• r > 0: exponential growth
• r < 0: exponential decay  
• r = 0: constant population

Applications:
• Early stages of bacterial growth
• Human population growth (historical periods)
• Radioactive decay (r < 0)
• Compound interest

Example: Bacteria doubling every 20 minutes
Doubling time: t_d = ln(2)/r
20 min = ln(2)/r → r = ln(2)/20 ≈ 0.0347 min⁻¹

After 2 hours: P(120) = P₀e^(0.0347×120) = P₀e^4.16 ≈ 64P₀
Population increases 64-fold

Logistic Growth Model

Assumptions:
• Limited carrying capacity
• Competition for resources increases with population
• Growth rate decreases as population approaches limit

Model: dP/dt = rP(1 - P/K)

Where K = carrying capacity (maximum sustainable population)

Solution: P(t) = K/(1 + Ae^(-rt)) where A = (K - P₀)/P₀

Properties:
• S-shaped (sigmoid) growth curve
• Inflection point at P = K/2
• Approaches K asymptotically
• Maximum growth rate occurs at P = K/2

Phase Analysis:
• P < K: dP/dt > 0 (population increases)
• P = K: dP/dt = 0 (equilibrium)
• P > K: dP/dt < 0 (population decreases)

Applications:
• Population growth with resource constraints
• Technology adoption (S-curves)
• Tumor growth
• Market penetration

Example: Yeast Population
r = 0.5 day⁻¹, K = 1000, P₀ = 10

P(t) = 1000/(1 + 99e^(-0.5t))

At t = 0: P = 10
At t = 5: P ≈ 76
At t = 10: P ≈ 607  
At t = 15: P ≈ 950
As t → ∞: P → 1000

Predator-Prey Models

Lotka-Volterra Model

System of two species: prey (x) and predator (y)

Assumptions:
• Prey grows exponentially without predators
• Predation rate proportional to prey × predator
• Predator deaths proportional to predator population
• Conversion efficiency from prey to predators constant

Model:
dx/dt = ax - bxy  (prey equation)
dy/dt = cxy - dy  (predator equation)

Where:
• a = prey growth rate
• b = predation efficiency  
• c = conversion efficiency
• d = predator death rate

Equilibrium: x* = d/c, y* = a/b

Properties:
• Closed orbits around equilibrium (periodic solutions)
• Conserved quantity: H = cx - d ln(x) + by - a ln(y)
• Oscillations with period ≈ 2π/√(ad)

Example: Lynx-Snowshoe Hare
Historical data shows ~10-year cycles
Parameters fitted to Canadian fur trade data

Modified Predator-Prey Models

Logistic Prey Growth:
dx/dt = ax(1 - x/K) - bxy
dy/dt = cxy - dy

Adds carrying capacity K for prey

Functional Response:
dx/dt = ax - (bxy)/(1 + bx/h) 
dy/dt = c(bxy)/(1 + bx/h) - dy

Saturating predation rate (Holling Type II)
h = handling time per prey item

Alternative Prey:
dx/dt = ax - (bxy)/(x + x₀)
dy/dt = c(bxy)/(x + x₀) - dy

Predators switch to alternative when main prey scarce

Age-Structured Models

Leslie Matrix Model for discrete age classes:

n(t+1) = L n(t)

Where n(t) is vector of population in each age class
L is Leslie matrix:

L = [f₁  f₂  f₃  ...  ]
    [s₁  0   0   ...  ]
    [0   s₂  0   ...  ]
    [⋮   ⋮   ⋱   ⋮   ]

fᵢ = fertility of age class i
sᵢ = survival probability from class i to i+1

Long-term Growth Rate:
λ = dominant eigenvalue of L
• λ > 1: population grows
• λ = 1: stable population  
• λ < 1: population declines

Stable Age Distribution:
Given by dominant eigenvector of L

Example: Human Population
Age classes: 0-14, 15-44, 45-64, 65+
Fertility mainly in 15-44 group
Different survival rates by age

Epidemic Models

SIR Model

Compartmental model dividing population into:
• S(t) = Susceptible individuals
• I(t) = Infected (and infectious) individuals  
• R(t) = Recovered (immune) individuals

Assumptions:
• Total population N = S + I + R is constant
• Infection rate proportional to S × I (mass action)
• Recovery rate proportional to I
• No births, deaths, or migration during epidemic

Model:
dS/dt = -βSI/N
dI/dt = βSI/N - γI  
dR/dt = γI

Where:
• β = transmission rate
• γ = recovery rate
• R₀ = β/γ = basic reproduction number

Key Insights:
• Epidemic occurs if R₀ > 1
• Peak infection: I_max occurs when S = N/R₀
• Final size relation: S(∞) = S₀e^(-R₀(1-S(∞)/N))

Example: COVID-19 Early Dynamics
β ≈ 0.5 day⁻¹, γ ≈ 0.1 day⁻¹
R₀ ≈ 5, epidemic spreads rapidly
Interventions reduce β (masks, distancing)

SEIR Model

Adds Exposed (E) compartment for incubation period:

dS/dt = -βSI/N
dE/dt = βSI/N - σE
dI/dt = σE - γI
dR/dt = γI

Where σ = 1/(incubation period)

More realistic for diseases with significant incubation time

SIS Model

For diseases without immunity (S → I → S):

dS/dt = γI - βSI/N
dI/dt = βSI/N - γI

With S + I = N, reduces to single equation:
dI/dt = I(β - γ) - βI²/N

Equilibrium: I* = N(1 - 1/R₀) if R₀ > 1, else I* = 0

Applications:
• Common cold
• Bacterial infections  
• Some sexually transmitted diseases

Network Models

Account for heterogeneous mixing patterns:

Rather than uniform mixing (mass action), consider:
• Contact networks (who interacts with whom)
• Spatial structure  
• Risk groups with different transmission rates

Example: Core Group Model
High-risk core group with high transmission
Low-risk population with lower transmission
Disease maintained in core, spills over to general population

Network Effects:
• Clustering: friends of friends are friends
• Degree distribution: some individuals have many contacts
• Small world: short path lengths despite clustering

Mechanical and Physical Models

Oscillatory Systems

Simple Harmonic Oscillator

Spring-mass system: F = ma
m(d²x/dt²) = -kx

Solution: x(t) = A cos(ωt + φ)
Where ω = √(k/m), A and φ from initial conditions

Properties:
• Period T = 2π/ω independent of amplitude
• Energy conserved: E = ½kA²
• Circular motion in phase space (x, dx/dt)

Applications:
• Mechanical springs
• Electrical LC circuits  
• Molecular vibrations
• Pendulum (small angles)

Damped Oscillator

Add damping: m(d²x/dt²) + c(dx/dt) + kx = 0

Characteristic equation: mr² + cr + k = 0
Discriminant Δ = c² - 4mk

Three cases:
1. Underdamped (Δ < 0): Oscillatory decay
   x(t) = Ae^(-γt) cos(ω_d t + φ)
   where γ = c/(2m), ω_d = √(ω₀² - γ²)

2. Critically damped (Δ = 0): Fastest return to equilibrium
   x(t) = (A + Bt)e^(-γt)

3. Overdamped (Δ > 0): Exponential decay without oscillation
   x(t) = Ae^(r₁t) + Be^(r₂t) where r₁, r₂ < 0

Quality factor: Q = ω₀/(2γ)
• High Q: lightly damped, many oscillations
• Low Q: heavily damped, few oscillations

Applications:
• Car suspensions (want critical damping)
• Building earthquake response
• Electronic circuits

Forced Oscillator

External driving force: m(d²x/dt²) + c(dx/dt) + kx = F₀cos(Ωt)

Steady-state solution: x(t) = A cos(Ωt - φ)

Amplitude: A = F₀/m / √[(ω₀² - Ω²)² + (2γΩ)²]
Phase lag: tan(φ) = 2γΩ/(ω₀² - Ω²)

Resonance occurs when Ω ≈ ω₀:
• Maximum amplitude transfer
• 90° phase lag at resonance
• Amplitude inversely proportional to damping

Applications:
• Musical instruments
• Tacoma Narrows Bridge collapse
• Atomic force microscopy
• Magnetic resonance imaging

Heat Transfer Models

Heat Equation

One-dimensional heat conduction:
∂u/∂t = α ∂²u/∂x²

Where:
• u(x,t) = temperature
• α = thermal diffusivity = k/(ρc)
• k = thermal conductivity
• ρ = density  
• c = specific heat

Boundary conditions:
• Prescribed temperature: u(0,t) = T₁, u(L,t) = T₂
• Insulated end: ∂u/∂x = 0
• Convective cooling: ∂u/∂x = h(u - T_ambient)

Solution methods:
• Separation of variables
• Fourier series
• Green's functions
• Numerical methods

Example: Rod with insulated ends
Initial condition: u(x,0) = f(x)
Solution: u(x,t) = Σ A_n cos(nπx/L) exp(-α(nπ/L)²t)

Heat diffuses, temperature profile smooths exponentially

Newton’s Law of Cooling

Lumped parameter model for cooling object:
dT/dt = -h(T - T_ambient)

Solution: T(t) = T_ambient + (T₀ - T_ambient)e^(-ht)

Where h depends on:
• Surface area
• Heat transfer coefficient  
• Mass and specific heat

Time constant τ = 1/h characterizes cooling rate

Applications:
• Cooling coffee
• Body temperature after death (forensics)
• Electronic component cooling
• Food safety (refrigeration)

Example: Hot coffee cooling
Initial temperature: 85°C
Room temperature: 20°C  
h = 0.1 min⁻¹

After 10 minutes: T = 20 + 65e^(-1) ≈ 44°C
After 30 minutes: T ≈ 23°C

Chemical Reaction Kinetics

First-Order Reactions

A → B with rate proportional to [A]:
d[A]/dt = -k[A]

Solution: [A](t) = [A]₀e^(-kt)

Half-life: t₁/₂ = ln(2)/k (independent of initial concentration)

Example: Radioactive decay
¹⁴C → ¹⁴N + e⁻ + ν̄
Half-life ≈ 5730 years
Used for carbon dating artifacts

Second-Order Reactions

A + B → C with rate proportional to [A][B]:
d[A]/dt = -k[A][B]

If [A]₀ = [B]₀ = a:
[A](t) = a/(1 + akt)

If [A]₀ ≠ [B]₀:
More complex solution involving difference in initial concentrations

No constant half-life (depends on concentration)

Enzyme Kinetics (Michaelis-Menten)

E + S ⇌ ES → E + P

Assuming fast pre-equilibrium:
v = V_max[S]/(K_m + [S])

Where:
• v = reaction velocity
• V_max = maximum velocity (when enzyme saturated)
• K_m = Michaelis constant (substrate concentration at v = V_max/2)

Behavior:
• Low [S]: v ≈ (V_max/K_m)[S] (first-order)
• High [S]: v ≈ V_max (zero-order, saturated)

Lineweaver-Burk plot: 1/v vs 1/[S] gives straight line
Used to determine K_m and V_max experimentally

Economic Models

Supply and Demand

Basic Market Model

Demand: Q_d = a - bp  (decreases with price)
Supply: Q_s = c + dp  (increases with price)

Equilibrium: Q_d = Q_s
a - bp* = c + dp*
p* = (a - c)/(b + d)
Q* = (ad + bc)/(b + d)

Consumer surplus = ½b(p_max - p*)²
Producer surplus = ½d(p* - p_min)²

Effects of shifts:
• Demand increase (a ↑): both price and quantity increase
• Supply increase (c ↑): price decreases, quantity increases
• Tax on suppliers: supply curve shifts up by tax amount

Market Dynamics

Price adjustment mechanism:
dp/dt = α(Q_d(p) - Q_s(p))

Where α > 0 is speed of price adjustment

Stability: If demand slope < 0 and supply slope > 0,
equilibrium is stable

Cobweb model (discrete time):
p_{t+1} based on Q_s(p_t)
Q_{t+1} based on p_{t+1}

Can lead to:
• Convergence to equilibrium
• Oscillations
• Divergent behavior

Growth Models

Solow Growth Model

Production function: Y = F(K, L) = AK^α L^(1-α)
Where:
• Y = output
• K = capital stock  
• L = labor
• A = technology parameter
• α = capital's share (0 < α < 1)

Capital accumulation: dK/dt = sY - δK
Where s = savings rate, δ = depreciation rate

In per capita terms (k = K/L, y = Y/L):
dk/dt = sy - (n + δ)k

Steady state: k* = [sA/(n + δ)]^(1/(1-α))

Properties:
• Countries converge to steady-state income
• Growth rate depends on technology growth
• Savings rate affects level, not long-term growth rate

Innovation and Endogenous Growth

AK model: Y = AK (no diminishing returns to capital)
Growth rate: g = sA - δ

R&D-driven growth:
dA/dt = λA(research effort)

Human capital models:
Include education and skill development
Growth can be sustained through human capital accumulation

Network effects and increasing returns:
Value increases with number of users
Can lead to multiple equilibria and path dependence

Financial Models

Black-Scholes Option Pricing

Stock price follows geometric Brownian motion:
dS = μS dt + σS dW

Where:
• μ = drift rate
• σ = volatility
• dW = random increment (Wiener process)

Option price V(S,t) satisfies:
∂V/∂t + ½σ²S² ∂²V/∂S² + rS ∂V/∂S - rV = 0

Boundary conditions for call option:
V(S,T) = max(S - K, 0)  at expiration
V(0,t) = 0  if stock worthless
V(S,t) ≈ S  for very large S

Solution: V = SN(d₁) - Ke^(-r(T-t))N(d₂)

Where N(x) is cumulative standard normal distribution
d₁ = [ln(S/K) + (r + σ²/2)(T-t)] / [σ√(T-t)]
d₂ = d₁ - σ√(T-t)

Greeks measure sensitivities:
• Delta: ∂V/∂S (hedge ratio)
• Gamma: ∂²V/∂S²
• Theta: ∂V/∂t (time decay)
• Vega: ∂V/∂σ

Portfolio Optimization

Markowitz mean-variance optimization:
Minimize risk: ½w^T Σw
Subject to: expected return ≥ target
            Σw_i = 1 (weights sum to 1)

Where:
• w = vector of portfolio weights
• Σ = covariance matrix of returns
• μ = vector of expected returns

Efficient frontier: trade-off between risk and return
Capital Asset Pricing Model (CAPM):
Expected return = risk-free rate + β × market risk premium

Risk parity:
Each asset contributes equally to portfolio risk
w_i × (Σw)_i = constant for all i

Biological Models

Pharmacokinetics

Single Compartment Model

Drug absorption and elimination:
dA/dt = -k_e A  (first-order elimination)

For oral dose D with absorption:
A(t) = (Dk_a/(k_a - k_e))[e^(-k_e t) - e^(-k_a t)]

Where:
• k_a = absorption rate constant
• k_e = elimination rate constant

Plasma concentration: C(t) = A(t)/V_d
Where V_d = apparent volume of distribution

Key parameters:
• Half-life: t₁/₂ = ln(2)/k_e
• Clearance: Cl = k_e V_d
• Bioavailability: F = fraction absorbed

Applications:
• Dosing regimens
• Drug interactions
• Individual patient optimization

Multi-Compartment Models

Two-compartment model:
Central compartment (blood) ⇌ Peripheral compartment (tissues)

dA₁/dt = -(k₁₀ + k₁₂)A₁ + k₂₁A₂
dA₂/dt = k₁₂A₁ - k₂₁A₂

Solution shows bi-exponential decay:
C(t) = Ae^(-αt) + Be^(-βt)

Where α > β are disposition rate constants
α represents distribution phase (fast)
β represents elimination phase (slow)

Used for:
• Antibiotics
• Anesthetics  
• Chemotherapy agents

Allometric Scaling

Biological quantities scale with body mass:
Y = aM^b

Where:
• Y = biological quantity (metabolic rate, lifespan, etc.)
• M = body mass
• a = proportionality constant
• b = scaling exponent

Examples:
• Basal metabolic rate: b ≈ 0.75
• Heart rate: b ≈ -0.25  
• Lifespan: b ≈ 0.25
• Brain mass: b ≈ 0.75

Kleiber's law: Metabolic rate ∝ M^(3/4)
Surface area ∝ M^(2/3), so metabolism not limited by surface area

Applications:
• Drug dosing across species
• Ecological energetics
• Evolutionary biology
• Biomedical scaling

Neural Models

Hodgkin-Huxley Model

Action potential in neurons:
C dV/dt = I - I_Na - I_K - I_L

Where:
• V = membrane potential
• I = applied current
• I_Na, I_K, I_L = ionic currents

Sodium current: I_Na = g_Na m³h(V - E_Na)
Potassium current: I_K = g_K n⁴(V - E_K)
Leak current: I_L = g_L(V - E_L)

Gating variables m, h, n follow:
dm/dt = α_m(V)(1-m) - β_m(V)m

Rate functions α, β are voltage-dependent
Derived from voltage clamp experiments

Model explains:
• Action potential shape and timing
• Refractory period
• Threshold behavior
• Propagation velocity

Simplified models (FitzHugh-Nagumo, Integrate-and-fire)
capture essential dynamics with fewer variables

Environmental and Climate Models

Greenhouse Effect

Energy balance model:
Incoming solar radiation = Outgoing thermal radiation

(1 - α)S₀/4 = εσT⁴

Where:
• α = planetary albedo (reflectivity)
• S₀ = solar constant ≈ 1367 W/m²
• ε = effective emissivity
• σ = Stefan-Boltzmann constant
• T = effective temperature

Without atmosphere: T ≈ 255 K (-18°C)
With greenhouse gases: T ≈ 288 K (15°C)
Greenhouse effect: ΔT ≈ 33 K

Enhanced greenhouse effect from CO₂:
ΔT = λ Δln(CO₂)
Where λ ≈ 3 K per CO₂ doubling (including feedbacks)

Feedbacks:
• Water vapor (positive): warmer air holds more moisture
• Ice-albedo (positive): less ice → less reflection → more warming
• Cloud feedback (uncertain): complex effects on radiation balance

Carbon Cycle

Box model with reservoirs:
• Atmosphere: ~850 GtC
• Ocean surface: ~1000 GtC  
• Deep ocean: ~38,000 GtC
• Vegetation: ~550 GtC
• Soils: ~1600 GtC

Flux equations:
dC_atm/dt = Emissions + Ocean→Atm - Atm→Ocean + Land→Atm - Atm→Land

Ocean-atmosphere exchange:
F = k(C_atm - αC_ocean)
Where α accounts for solubility and chemistry

Land-atmosphere exchange:
Photosynthesis: F_photo = β(CO₂)G(climate)
Respiration: F_resp = Q₁₀^((T-T₀)/10) × R₀

Where Q₁₀ ≈ 2 is temperature sensitivity of respiration

Airborne fraction: ~45% of emissions remain in atmosphere
Ocean and land act as carbon sinks, but saturation possible

Population Ecology

Metapopulation Models

Population in patches connected by migration:
dp_i/dt = c_i p_i(1 - p_i) - e_i p_i + Σ_j m_ij p_j

Where:
• p_i = probability patch i is occupied
• c_i = colonization rate
• e_i = extinction rate
• m_ij = migration from patch j to i

Mainland-island model:
dp/dt = c(1 - p) - ep

Equilibrium: p* = c/(c + e)

Species persists if colonization > extinction

Spatial structure affects persistence:
• Corridor effects
• Edge effects
• Fragmentation impacts
• Conservation reserve design

Food Web Dynamics

Multi-species interactions:
dx_i/dt = f_i(x₁, x₂, ..., x_n)

Interaction matrix A with elements:
• a_ij > 0: species j benefits species i
• a_ij < 0: species j harms species i  
• a_ij = 0: no interaction

Stability analysis:
System stable if all eigenvalues of Jacobian have negative real parts

May's complexity-stability paradox:
Random large food webs typically unstable
But real food webs show structure that promotes stability:
• Weak interactions
• Compartmentalization  
• Non-random topology

Model Validation and Analysis

Sensitivity Analysis

How do model outputs change with parameter variations?

Local Sensitivity:
S_ij = ∂Y_i/∂p_j × p_j/Y_i (normalized partial derivatives)

Global Sensitivity:
Vary parameters over their full ranges
Use sampling methods (Monte Carlo, Latin Hypercube)

Sobol indices:
First-order: S_i = V(E(Y|X_i))/V(Y)
Total: S_Ti = E(V(Y|X_{~i}))/V(Y)

Where V = variance, E = expectation

Applications:
• Parameter estimation priorities
• Model simplification  
• Uncertainty quantification
• Robustness assessment

Example: Epidemiological model
If S_R₀ >> S_other, focus on estimating transmission parameters
If S_initial conditions high, need good initial data

Uncertainty Quantification

Sources of uncertainty:
• Parameter uncertainty
• Model structure uncertainty
• Input uncertainty
• Numerical uncertainty

Probability distributions for parameters:
p(θ) based on prior knowledge and data

Bayesian updating:
p(θ|data) ∝ p(data|θ) × p(θ)

Propagation methods:
• Monte Carlo sampling
• Polynomial chaos expansion  
• Ensemble methods

Confidence/credible intervals for predictions
Prediction intervals accounting for all uncertainties

Example: Climate projections
Multiple models, parameter ensembles
Structural uncertainty from different physics representations

Model Selection and Comparison

Information criteria for model comparison:
AIC = -2 log(likelihood) + 2k
BIC = -2 log(likelihood) + k log(n)

Where k = number of parameters, n = sample size

Lower values indicate better models
BIC penalizes complexity more heavily

Cross-validation:
Split data into training and testing sets
Fit model on training, evaluate on testing
k-fold CV: repeat with different splits

Bayes factors:
B₁₂ = p(data|M₁)/p(data|M₂)
Ratio of marginal likelihoods

Model averaging:
Weight predictions by model probabilities
Accounts for model uncertainty

Practical considerations:
• Purpose of model (prediction vs. understanding)
• Computational cost
• Interpretability
• Extrapolation requirements

Identifiability Analysis

Can model parameters be uniquely determined from data?

Structural identifiability:
Mathematical question - given perfect, noise-free data
Analyze using:
• Taylor series coefficients
• Similarity transformations
• Differential algebra

Practical identifiability:  
Given realistic data with noise and limited observations
Use:
• Fisher Information Matrix
• Profile likelihood
• Bootstrap confidence intervals

Non-identifiable parameters:
• Remove from model (structural)
• Fix at reasonable values
• Use prior information (Bayesian)
• Collect more informative data

Example: Two-parameter exponential model
y = ae^(-bx)
If only late-time data available, a and b not separately identifiable
Product ab may be identifiable

Computational Implementation

Model Development Workflow

1. Conceptual Model
   • Diagram relationships
   • State assumptions clearly
   • Define variables and parameters

2. Mathematical Formulation  
   • Write equations
   • Specify initial/boundary conditions
   • Identify dimensionless groups

3. Analytical Analysis
   • Equilibria and stability
   • Limiting cases
   • Scaling behavior

4. Numerical Implementation
   • Choose appropriate algorithms
   • Implement in software
   • Verify implementation

5. Validation and Testing
   • Compare with analytical solutions
   • Test limiting cases
   • Sensitivity analysis

6. Application and Interpretation
   • Run simulations
   • Interpret results
   • Communicate findings

Software Tools

General Purpose:
• MATLAB: Extensive toolboxes, excellent for prototyping
• Python: NumPy/SciPy for numerics, Matplotlib for plotting
• R: Statistical analysis and modeling
• Mathematica: Symbolic and numerical computation

Specialized:
• Systems Biology Markup Language (SBML)
• COPASI: Biochemical network simulation
• NetLogo: Agent-based modeling
• AnyLogic: Multi-method modeling platform

High Performance:
• Julia: High-performance scientific computing
• C/C++/Fortran: Maximum speed, more development time
• GPU computing: CUDA, OpenCL for parallel problems

Version Control and Reproducibility:
• Git for version control
• Jupyter notebooks for documentation
• Container systems (Docker) for reproducible environments

Best Practices

Code Organization:
• Separate model definition from analysis
• Use meaningful variable names
• Comment assumptions and derivations
• Modular design for reusability

Documentation:
• Model purpose and assumptions
• Parameter meanings and units  
• Validation tests performed
• Known limitations

Reproducibility:
• Fixed random seeds
• Version control
• Dependency management
• Clear instructions for running

Quality Assurance:
• Unit tests for key functions
• Regression tests
• Peer review
• Open source when possible

Example Python structure:
model/
├── core.py          # Model equations
├── parameters.py    # Parameter definitions  
├── analysis.py      # Analysis functions
├── validation.py    # Tests and validation
├── examples/        # Example applications
└── docs/           # Documentation

Case Studies

COVID-19 Modeling

Challenge: Rapid response to emerging pandemic

Initial Models (March 2020):
• Simple SIR with estimated R₀ ≈ 2.5
• Exponential growth phase analysis
• Crude estimates of hospital capacity needs

Model Evolution:
• SEIR with latent period
• Age-structured models (different risks)  
• Spatial models (geographic spread)
• Network models (superspreading)

Key Insights:
• Importance of asymptomatic transmission
• Role of superspreading events
• Impact of interventions (lockdowns, masks)
• Waning immunity and variants

Challenges:
• Rapidly changing virus characteristics
• Behavioral responses to policies
• Data quality and reporting delays
• Communication of uncertainty

Lessons Learned:
• Multiple models better than single model
• Importance of real-time data
• Need for scenario planning
• Value of transparent communication

Climate Change Projections

Problem: Long-term climate prediction with multiple uncertainties

Model Hierarchy:
• Energy Balance Models: Global temperature
• Box Models: Carbon cycle, ocean circulation
• General Circulation Models: Full 3D dynamics

Key Components:
• Radiative forcing from greenhouse gases
• Climate sensitivity (temperature response)
• Carbon cycle feedbacks
• Ice sheet dynamics

Major Uncertainties:
• Cloud feedbacks
• Carbon cycle responses  
• Regional precipitation patterns
• Tipping points and nonlinearities

IPCC Approach:
• Multiple model ensembles
• Scenario-based projections (RCP, SSP)
• Probabilistic statements
• Regular assessment reports

Policy Applications:
• Global temperature targets (1.5°C, 2°C)
• Carbon budgets for emissions
• Regional impact assessments
• Adaptation planning

Success Factors:
• International collaboration
• Model intercomparison projects
• Open data and code sharing
• Regular model evaluation

Drug Development

Application: Optimize dosing regimens for new medications

Pharmacokinetic Modeling:
• ADME: Absorption, Distribution, Metabolism, Excretion
• Population PK: Account for individual variability
• Physiologically-based PK: Mechanistic organ models

Pharmacodynamic Modeling:
• Dose-response relationships
• Efficacy and safety endpoints
• Time-course of drug effects

Model-Informed Drug Development:
• Optimal dose selection
• Clinical trial design
• Regulatory submissions

Example: Oncology Drug
Population: Cancer patients with varying organ function
Challenge: Narrow therapeutic window (efficacy vs. toxicity)

Model Components:
• 2-compartment PK with nonlinear elimination
• Tumor growth inhibition model
• Hematological toxicity model

Covariates:
• Body weight, age, organ function
• Disease severity
• Concomitant medications

Simulation Results:
• Dose recommendations by patient subgroup
• Probability of target attainment
• Risk-benefit analysis

Regulatory Acceptance:
• FDA/EMA guidance on M&S
• Model qualification process
• Integration with clinical data

Practice Problems

Population Dynamics

1. A bacterial culture starts with 100 cells and doubles every 30 minutes.
   a) Write the exponential growth model
   b) How many cells after 3 hours?
   c) When will population reach 10,000?

2. A population follows logistic growth with r = 0.1/year and K = 1000.
   If current population is 200:
   a) What is the growth rate now?
   b) When will population reach 500?
   c) What is maximum growth rate and when does it occur?

3. Predator-prey system: dx/dt = x - xy, dy/dt = -y + xy
   a) Find equilibrium point
   b) Classify stability
   c) Find conserved quantity

Physical Models

1. A 95°C cup of coffee cools according to Newton's law with h = 0.05 min⁻¹.
   Room temperature is 20°C.
   a) Write the differential equation
   b) What is temperature after 10 minutes?
   c) How long until temperature is 40°C?

2. Damped oscillator: d²x/dt² + 2dx/dt + 5x = 0
   a) Classify the damping
   b) Find the solution with x(0) = 1, x'(0) = 0
   c) What is the period of oscillation?

3. Heat equation on [0,1] with u(0,t) = u(1,t) = 0, u(x,0) = sin(πx)
   a) Use separation of variables
   b) Write the solution
   c) What is temperature at x = 0.5 after time t = 0.1?

Economics/Biology

1. Supply: Qs = 2p, Demand: Qd = 10 - p
   a) Find equilibrium price and quantity
   b) What happens if government sets price at $2?
   c) Calculate consumer and producer surplus at equilibrium

2. SIR epidemic model with β = 0.3, γ = 0.1, N = 1000
   Initial conditions: S(0) = 990, I(0) = 10, R(0) = 0
   a) What is R₀?
   b) Will epidemic occur?  
   c) What is final size of epidemic?

3. Drug concentration follows C(t) = 100e⁻⁰·⁵ᵗ mg/L
   a) What is elimination half-life?
   b) What is concentration after 4 hours?
   c) When does concentration drop below 10 mg/L?

Answer Key

Population Dynamics

1. a) P(t) = 100 × 2^(t/0.5) = 100 × 2^(2t) where t in hours
   b) P(3) = 100 × 2^6 = 6,400 cells  
   c) 10,000 = 100 × 2^(2t) → t = log₂(100)/2 ≈ 3.32 hours

2. a) dP/dt = 0.1 × 200 × (1 - 200/1000) = 16/year
   b) Solve 500 = 1000/(1 + Ae^(-0.1t)), get t ≈ 10.99 years
   c) Max at P = K/2 = 500, rate = rK/4 = 25/year

3. a) Equilibrium: (1, 1)
   b) Jacobian eigenvalues: ±i, center (neutrally stable)
   c) H = x - ln(x) + y - ln(y) = constant

Physical Models

1. a) dT/dt = -0.05(T - 20), T(0) = 95
   b) T(10) = 20 + 75e^(-0.5) ≈ 65.5°C
   c) 40 = 20 + 75e^(-0.05t) → t = 20 ln(3.75) ≈ 26.1 min

2. a) Discriminant = 4 - 20 = -16 < 0, underdamped
   b) x(t) = e^(-t)[cos(2t) + 0.5sin(2t)]
   c) ωd = 2, period = π seconds

3. a) u(x,t) = Σ An sin(nπx)e^(-n²π²t)
   b) u(x,t) = sin(πx)e^(-π²t) (n=1 term only)  
   c) u(0.5, 0.1) = sin(π/2)e^(-π²×0.1) = e^(-π²/10) ≈ 0.00005

Economics/Biology

1. a) 2p = 10 - p → p* = 10/3, Q* = 20/3
   b) At p = $2: shortage of 6 units (Qd = 8, Qs = 4)
   c) CS = ½×(20/3)×(10-10/3) = 200/9, PS = ½×(20/3)×(10/3) = 100/9

2. a) R₀ = β/γ = 0.3/0.1 = 3
   b) Yes, since R₀ > 1
   c) Final size equation gives S(∞) ≈ 25, so ~975 infected total

3. a) Half-life = ln(2)/0.5 = 1.39 hours
   b) C(4) = 100e^(-2) ≈ 13.5 mg/L
   c) 10 = 100e^(-0.5t) → t = 2ln(10) ≈ 4.61 hours

Quick Reference

Model Types

Deterministic vs. Stochastic: Fixed vs. random components
Discrete vs. Continuous: Time/space discretization
Linear vs. Nonlinear: Superposition principle
Static vs. Dynamic: Equilibrium vs. time evolution

Key Techniques

Conservation Laws: Mass, energy, momentum conservation
Dimensional Analysis: Scaling and similarity
Equilibrium Analysis: Steady states and stability
Perturbation Methods: Small parameter expansions

Validation Methods

Sensitivity Analysis: Parameter importance
Uncertainty Quantification: Error propagation
Model Selection: Information criteria, cross-validation
Identifiability: Parameter estimability

Common Models

Population: Exponential, logistic, predator-prey
Epidemic: SIR, SEIR network models
Physical: Oscillators, heat equation, reaction kinetics
Economic: Supply/demand, growth models, option pricing