Modelli Matematici applicati all`Ecologia

S. Console - M. Roggero
Modelli Matematici
applicati all’Ecologia
LAUREA MAGISTRALE in:
ANALISI E GESTIONE DELL’AMBIENTE
(AGAm)
A.A. 2006/2007
Indice
Capitolo 1 - Introduzione ai modelli matematici
Modelli statici . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Modelli dinamici . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Capitolo 2 - Equazioni differenziali
Equazioni differenziali del primo ordine . . . . . . . . . . . . . .
Equazioni differenziali del primo ordine e campi di direzioni
Equazioni differenziali a variabili separabili . . . . . . . . .
Equazioni differenziali del secondo ordine . . . . . . . . . . . . .
Equazioni del 2° ordine omogenee a coefficienti costanti . .
Cenni sulla formula di Taylor . . . . . . . . . . . . . . . . . . . .
Esercizi risolti . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Altri esercizi . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Capitolo 3 - Modelli di popolazioni isolate
Modelli di popolazioni isolate omogenee . . . . . . . .
Modello Malthusiano . . . . . . . . . . . . . . . .
Il modello logistico . . . . . . . . . . . . . . . . .
Modelli di popolazioni divise in due classi d’età . . . .
Le matrici . . . . . . . . . . . . . . . . . . . . . . . . .
Operazioni tra matrici . . . . . . . . . . . . . . .
Il determinante di una matrice quadrata . . . . .
Riduzione di una matrice. . . . . . . . . . . . . .
Sistemi lineari . . . . . . . . . . . . . . . . . . . .
Autovalori e autovettori di una matrice quadrata
Potenze di una matrice . . . . . . . . . . . . . . .
Modelli di popolazioni strutturate per età o per taglia
Esercizi . . . . . . . . . . . . . . . . . . . . . . . . . .
Spunti per approfondimenti . . . . . . . . . . . . . . .
Altri modelli per popolazioni isolate . . . . . . .
Tabelle e grafi di vita e matrici di Leslie . . . . .
2
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
5
6
6
.
.
.
.
.
.
.
.
8
8
9
10
11
12
13
14
14
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
16
17
18
19
24
26
27
30
31
33
34
34
35
36
38
38
40
Indice
3
Capitolo 4 - Calibrazione di modelli matematici
Trasformazioni linearizzanti . . . . . . . . . . . . . . . . . . . . . . . . . .
Calibrazione per un modello di crescita logistica . . . . . . . . . . . . . .
Esercizi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
44
48
48
53
Capitolo 5 - Sistemi di equazioni differenziali
Il piano delle fasi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Soluzioni di sistemi lineari omogenei a coefficienti costanti . . . . . . . . .
Esercizi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
55
57
60
64
Capitolo 6 - Stabilità
Punti di equilibrio di un’equazione differenziale . . . . . . . . . . . .
Equilibrio: Stabile o instabile? . . . . . . . . . . . . . . . . . .
Stabilità per i punti di equilibrio di una equazione differenziale
Stabilità e biforcazione . . . . . . . . . . . . . . . . . . . . . . .
Stabilità per gli equilibri di modelli discreti . . . . . . . . . . . . . .
Stabilità per punti di equilibrio di sistemi di equazioni differenziali .
Piano delle fasi e stabilità per sistemi lineari omogenei . . . . .
Stabilità per sistemi generali . . . . . . . . . . . . . . . . . . .
Esercizi risolti . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Altri esercizi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
65
65
66
66
67
70
76
76
82
83
84
Capitolo 7 - Modelli di popolazioni conviventi
Modelli quadratici . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Predazione: il modello Lotka-Volterra . . . . . . . . . . . . . . . . . . . .
Studio delle equazioni di Lotka-Volterra . . . . . . . . . . . . . . . .
Un’applicazione del modello: effetto pesca . . . . . . . . . . . . . . .
Autolimitazione della preda . . . . . . . . . . . . . . . . . . . . . . .
Mutua interferenza . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Risposta funzionale del predatore . . . . . . . . . . . . . . . . . . . .
Competizione interspecifica . . . . . . . . . . . . . . . . . . . . . . . . . .
Modello Lotka-Volterra quadratico per la competizione di due specie
Cooperazione o mutualismo . . . . . . . . . . . . . . . . . . . . . . . . . .
Modello quadratico per la cooperazione obbligatoria . . . . . . . . .
Modello quadratico per la cooperazione facoltativa . . . . . . . . . .
Limite dell’applicazione del modello quadratico per il mutualismo . .
85
87
89
90
95
98
98
98
101
101
104
108
109
111
Appendice 1 - Le successioni
Successioni aritmetiche . . . .
Successioni geometriche . . .
Applicazioni delle successioni
Esercizi risolti . . . . . . . . .
116
117
118
119
121
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
Modelli Matematici applicati all’Ecologia (28.07.06)
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
4
Indice
Altri esercizi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123
Appendice 2 - I numeri complessi
Definizione di C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Il teorema fondamentale dell’algebra . . . . . . . . . . . . . . . . . . . . .
Esercizi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
125
125
127
127
Appendice 3 - Le funzioni di due o più variabili
129
Funzioni a più variabili e loro rappresentazione grafica . . . . . . . . . . . 129
Derivate parziali . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130
S. Console – M. Roggero
Capitolo 1
Introduzione ai modelli matematici
Un modello matematico è una rappresentazione astratta di un sistema, cioè un
insieme di processi evolutivi interconnessi.
La forma del modello dipende dagli scopi per cui viene costruito e dalle conoscenze disponibili (dati sperimentali e conoscenza del suo funzionamento).
Molto schematicamente possiamo dire che ci sono due tipi di modelli: deterministici e statistici.
In genere si sceglie un modello statistico se non si conoscono le leggi che governano
il funzionamento del sistema o se la complessità del sistema è tale per cui non è (anche
computazionalmente) possibile applicare uno schema deterministico.
Tuttavia il confine tra
i due modelli non è cosı̀
netto. Si pensi al semplice esempio di un tavolo da biliardo. Se sul tavolo c’è una sola pallina,
allora il moto della palla
può essere descritto secondo un modello deterministico: se conosco la velocità e la direzione in cui la colpisco, semplici leggi fisiche
mi permettono di prevederne la traiettoria e quindi di costruire un modello deterministico. Se invece sul tavolo ho molte palline, allora il sistema è cosı̀ complesso da
non poter essere descritto da un modello deterministico.
In questo corso ci occuperemo solo di modelli deterministici.
All’interno di tali modelli possiamo fare un’ulteriore distinzione tra modelli statici
e modelli dinamici.
5
6
§ 1.1
Capitolo 1
Modelli statici
Sono modelli che non dipendono dal tempo e sono quindi costituiti da equazioni algebriche. Un semplice esempio è il seguente.
Consideriamo due canali di portate Q1 e Q2
che confluiscono in una canale di portata Qv .
Supponiamo che i canali contengano una sostanza disciolta con concentazioni C1 , C2 e Cv .
Allora la portata Qv e la concentazione Cv
sono regolati da semplici equazioni algebriche.
Il bilancio delle portate porta all’equazione
Qv = Q1 + Q2 .
Il bilancio delle sostenze disciolte
Qv Cv = Q1 C1 + Q2 C2 ,
permette di ottenere
Cv =
§ 1.2
Q1 C1 + Q2 C2
.
Q1 + Q2
Modelli dinamici
Si tratta di modelli che descrivono l’evoluzione del sistema nel tempo. Nel caso
più semplice, lo stato del sistema al tempo t è descritto da una funzione x(t). Ad
esempio, in dinamica delle popolazioni, x(t) rappresenta la popolazione presente al
tempo t. L’equazione che regola il sistema è dinamica: si descrive la variazione
istantanea della popolazione (modello continuo) oppure la variazione in intervalli
fissi di tempo (modello tempo-discreto).
Esempio 1.1 (crescita esponenziale, modello continuo). Supponiamo che la
crescita istantanea di una certa popolazione (ad esempio una coltura di batteri) oppure il diffondersi di una malattia sia proporzionale alla popolazione presente. In
tal caso il modello è costituito dall’equazione differenziale
dx
= kx(t) ,
dt
con k tasso di crescita costante.
Le soluzioni dell’equazione differenziale sono date da
x(t) = x0 ekt ,
S. Console – M. Roggero
Introduzione ai modelli matematici
7
con x0 che rappresenta la popolazione presente al tempo t = 0. Tale modello è stato
per primo proposto dall’economista Malthus alla fine del 1700 per spiegare la crescita
di popolazione durante la rivoluzione industriale. In figura si mostra il grafico di x(t)
nel caso in cui k > 0 (a sinistra) e k < 0 (a destra).
Esempio 1.2 (crescita esponenziale, modello discreto). Supponiamo che una
popolazione cresca in intervalli di tempi ben delimitati. Questo avviene ad esempio
per animali che depongono le uova solo una volta all’anno. Indichiamo con xn la
popolazione dopo n intervalli di tempo. Supponiamo che il tasso di crescita sia
costante k. allora
xn+1 − xn = kxn ,
con k la crescita costante.
Otteniamo quindi la successione ricorrente
x0 = a ,
xn+1 = (1 + k)xn ,
con a che rappresenta la popolazione iniziale. Osserviamo che la successione è in
effetti una progressione geometrica
xn = a(1 + k)n .
Modelli Matematici applicati all’Ecologia (3.12.06)
Capitolo 2
Equazioni differenziali
I modelli matematici per lo studio di una popolazione isolata sono equazioni
differenziali. Premettiamo dunque allo studio dei modelli di popolazioni isolate una
breve introduzione alle equazioni differenziali.
§ 2.1
Equazioni differenziali del primo ordine
Esprimono un legame tra
ˆ una variabile indipendente x;
ˆ una variabile dipendente y (che sta ad indicare una funzione y = y(x));
ˆ una o più derivate della funzione y.
Esempi 2.1.
1. y 0 = f (x) è un’equazione differenziale una cui soluzione è una
primitiva di y. Ad esempio se y 0 = ex + x2 , una soluzione è y = ex + x3 /3 ed
una generica soluzione è y = ex + x3 /3 + c, con c una costante arbitraria.
2. y 0 = y ha come soluzione y(x) = ex ed una generica soluzione è y(x) = cex ,
con c una costante arbitraria. Infatti la derivata di y(x) = cex è proprio y
stesso.
Dagli esempi vediamo che le soluzioni delle equazioni differenziali non sono uniche
ma dipendono da una o, in generale, piu’ costanti arbitrarie.
Una equazione differenziale in cui compare solo la derivata prima della funzione
incognita si dice equazione differenziale del 1o ordine; se compare anche y 00 si dice
del 2 ordine, e cosı̀ via. Un’equazione differenziale del primo ordine si dice in forma
normale se è del tipo:
y 0 = F (x, y) .
Esempio 2.2. y 0 = y è un’equazione differenziale del primo ordine in forma normale.
Un’equazione in forma normale si dice autonoma se è della forma y 0 = F (y), cioè
la funzione F non dipende da x (che rappresenta il tempo, in molte delle applicazioni
che vedremo).
8
Equazioni differenziali
9
Equazioni differenziali del primo ordine e campi di direzioni
Un’equazione differenziale del primo ordine in forma normale y 0 = F (x, y) può essere
interpretata geometricamente come l’associare a ogni punto (x, y) del piano cartesiano una direzione (determinata da y 0 che è il coefficiente angolare della funzione
incognita y = y(x)).
In altre parole ogni soluzione
y = y(x) dovrà avere in (x, y) retta tangente di coefficiente angolare y 0 = F (x, y). Si determina
allora un campo di direzioni associando a ogni (x, y) la direzione
determinata dal coefficiente angolare y 0 = F (x, y).
Ad esempio l’equazione (autonoma, che descrive il modello logistico) y 0 = y(1 − y) determina il
campo di direzioni in figura.
Intuitivamente, se partiamo da
un punto P e ci muoviamo seguendo punto per punto la direzione data dal campo di direzioni, la nostra traiettoria è
una curva; si tratta proprio di una delle soluzioni, l’unica che passa per P .
Questo giustifica (intuitivamente) il seguente
1.5
1
y(t)
0.5
0
1
2
t
3
4
–0.5
–1
Teorema 2.3 (di esistenza e unicità della soluzione di un’equazione differenziale del 1° ordine in forma normale). Data un’equazione differenziale
y 0 = F (x, y) e un punto P (x0 , y0 ), sotto opportune condizioni di regolarità della
funzione F (x, y) esiste una e una sola funzione y = f (x) soluzione dell’equazione
differenziale tale che y0 = f (x0 ). (cioè tale che la curva y = f (x) passa per il punto
P (x0 , y0 ).)
Data un’equazione diffenziale y 0 = F (x, y) il problema di trovare la soluzione
y = y(x) tale che y0 = f (x0 ) si dice problema di Cauchy.
Osserviamo che il teorema di esistenza e unicità significa in particolare che il
campo di direzioni (e la condizione iniziale) determinano univocamente la soluzione.
Il nostro punto di vista sarà non tanto quello di scrivere esplicitamente le soluzioni,
quanto di capire il loro comportamento in base a proprietà della funzione F (x, y).
Questo è l’approccio che avremo ad esempio nel capitolo sulla stabilità.
Vediamo però ora un metodo di risoluzione delle equazioni differenziali, che si
può applicare sempre ad equazioni autonome.
Modelli Matematici applicati all’Ecologia (3.12.06)
10
Capitolo 2
1.5
1.5
1
1.5
1
y(t)
1
y(t)
0.5
0
y(t)
0.5
1
2
t
3
0.5
0
4
1
2
t
3
0
4
–0.5
–0.5
–0.5
–1
–1
–1
1
2
t
3
4
(a) Soluzione passante per un
(b) Soluzioni per varie con-
(c) Soluzioni e campo di di-
punto
dizioni iniziali
rezioni
Equazioni differenziali a variabili separabili
Un’equazione differenziale della forma
y 0 = f (x)g(y)
si dice a variabili separabili.
Esempio 2.4. L’equazione
√
y 0 = sin(x) y
è a variabili separabili.
dy
Per risolverla scrivo y 0 come
e tratto questo termine come una frazione.
dx
Moltiplico i due membri per dx e porto quindi tutto ciò che contiene la y dalla
stessa parte di dy e tutto ciò che contiene x da quella di dx, ottenendo:
dy
= f (x)dx.
g(y)
Nell’esempio scrivo
dy
√
= sin(x) y
dx
e poi
dy
√ = sin(x)dx
y
Infine, integrando ambo i membri trovo la soluzione generale.
S. Console – M. Roggero
Equazioni differenziali
11
Nell’esempio
Z
dy
√ =
y
Z
sin(x)dx
√
2 y = − cos(x) + c
y=
1
c
− cos(x) +
2
2
2
.
Esempio 2.5. Anche l’equazione del 1° ordine omogenea a coefficienti costanti
y 0 − ky = 0
è a variabili separabili. Portando ky a secondo membro si ottiene infatti:
y 0 = ky.
Risolvendo col metodo sopra illustrato si trovano le soluzioni:
y = c ekx .
Si noti che k è la soluzione dell’equazione T − k = 0 ottenuta mettendo T al posto
di y 0 e 1 al posto di y.
§ 2.2
Equazioni differenziali del secondo ordine
Una equazione diffenrenziale in cui compare anche y 00 oltre eventualmente a y 0 e a
y e a funzioni note della variabile x, ma non derivate della y di ordine superiore, si
dice del 2° ordine. Si dirà in forma normale se è del tipo:
y 00 = F (x, y, y 0 ).
Anche in questo caso vale un risultato di esistenza e unicità:
Teorema 2.6. Data un’equazione differenziale y 00 = F (x, y, y 0 ) e inoltre un punto
P (x0 , y0 ) e un numero reale m, sotto opportune condizioni di regolarità della funzione F (x, y, y 0 ) esiste una e una sola funzione y = f (x) soluzione dell’equazione
differenziale tale che y0 = f (x0 ) e m = f 0 (x0 ) (cioè tale che la curva y = f (x) passa
per il punto P (x0 , y0 ) e in quel punto la sua retta tangente ha coefficiente angolare
m.)
In questo caso l’equazione differenziale ha infinite soluzioni dipendenti da 2
parametri.
Modelli Matematici applicati all’Ecologia (3.12.06)
12
Capitolo 2
Equazioni del 2° ordine omogenee a coefficienti costanti
Una equazione differenziale della forma:
y 00 + αy 0 + β = 0
con α, β ∈ R si dice omogenea del 2° ordine a coefficienti costanti.
In analogia con quanto visto nell’Esempio 2.5 cerchiamo eventuali soluzioni di
tipo esponenziale y = ekx . Derivando due volte y 0 = kekx e y 00 = k 2 ekx e sostituendo
otteniamo:
ekx k 2 + αk + β = 0.
Se k è soluzione dell’equazione T 2 + αT + β = 0, la funzione y = ekx sarà soluzione
dell’equazione differnziale. Distinguiamo a seconda che l’equazione T 2 + αT + β = 0
abbia discriminante ∆ = α2 − 4β positivo, nullo oppure negativo
ˆ se ∆ > 0 ci solno due soluzioni reali distinte k1 e k2 . Le soluzioni dell’equazione
differenziale sono al variare di c1 , c2 ∈ R:
y = c1 ek1 x + c2 ek2 x
ˆ se ∆ = 0 c’è una sola soluzione reale (doppia) k = −α/2. Le soluzioni
dell’equazione differenziale sono al variare di c1 , c2 ∈ R:
y = ekx (c1 + c2 x)
√
ˆ se ∆ < 0 si ponga a = −α/2 e b = ( −∆)/2. Le soluzioni dell’equazione
differenziale sono al variare di c1 , c2 ∈ R:
y = eax (c1 cos(bx) + c2 sin(bx))
Per verificare la correttezza delle soluzioni anche nel secondo e terzo caso, basta
derivare e sostituire.
Consideriamo ora il caso più generale di una equazione differenziale a coefficienti
costanti non omogenea ossia del tipo:
y 00 + αy 0 + βy = f (x).
(1)
Possiamo osservare che se y1 (x) è una sua soluzione e y0 (x) è una soluzione dell’equazione omogenea associata (ottenuta sostituendo 0 al posto di f (x)), allora la loro
somma è ancora una sua soluzione. Potremo allora ottenere tutte le soluzioni della (1) sommando una (qualunque) sua soluzione particolare con ciascuna soluzione
dell’omogenea associata.
S. Console – M. Roggero
Equazioni differenziali
13
Esempio 2.7. Una soluzione di y 00 + 4y = x è la funzione y0 (x) = 41 x. Risolvendo
l’omogenea associata y 00 +4y = 0 e sommando a ciascuna soluzione la funzione y0 (x)
si trovano tutte le soluzioni di y 00 + 4y = x:
1
y(x) = x + c1 cos(2x) + c2 sin(2x).
4
La risoluzione di una equazione differenziale a coefficienti costanti procede quindi
in due passi. Il metodo prima illustrato permette di risolvere facilmente l’equazione
omogenea associata; purtroppo però non esiste un algoritmo che permetta anche di
trovare una soluzione particolare dell’equazione differenziale completa.
Un metodo che talvolta funziona è quello di cercare una soluzione che sia dello
stesso tipo di f (x). Un caso particolarmente importante è quello in cui f (x) è
polinomiale: in tal caso c’è sempre una soluzione che è una funzione polinomiale
dello stesso grado.
In generale quando non si riesce a determinare una qualche soluzione esatta, se
ne può cercare una approssimata di tipo polinomiale approssiamando f (x) con un
polinomio opportuno mediante lo sviluppo di Taylor.
§ 2.3
Cenni sulla formula di Taylor
Lo sviluppo di Taylor serve per approssimare una funzione localmente (vicino a un
punto) con un polinomio.
Data una funzione f (x) (derivabile fino all’ordine n) e un punto x0 del dominio
il polinomio di Taylor di ordine n è
1
1
Pn (x) = f (x0 ) + f 0 (x0 )(x − x0 ) + f 00 (x0 )(x − x0 )2 + f 000 (x0 )(x − x0 )3 + ....+
2!
3!
1
+ !f (n) (x0 )(x − x0 )n .
n!
Notiamo che il polinomio di Taylor di grado 1 in x0 altro non è che l’approssimazione
della funzione con la retta tangente nel punto x0 .
Sotto la condizione che la funzione abbia derivate continue fino all’ordine n, la
differenza tra f (x) e Pn (x) e’ trascurabile rispetto a (x − x0 )n (resto della formula
di Taylor secondo Peano). Intuitivamente, al crescere di n, il polinomio Pn (x) è
un’approssimazione sempre migliore di f (x) vicino a x0 .
Esempio 2.8. Consideriamo f (x) = sin(x). I suoi polinomi di Taylor calcolati nel
x3
punto x0 = 0 sono: P1 (x) = P2 (x) = x, P3 (x) = P4 (x) = x − , P5 (x) = P6 (x) =
6
x3 x5
x3 x5 x7
x−
+ , P7 (x) = P8 (x) = x −
+
− , ...
6
5!
6
5!
7!
Modelli Matematici applicati all’Ecologia (3.12.06)
14
Capitolo 2
Vediamo i disegni della funzione sin(x) e di P1 (x) (a
sinistra) e della funzione assieme ai polinomi di Taylor in 0, P1 (x) e P3 (x) (a
destra).
Se disegnamo ora la funzione sin(x) assieme ai polinomi di Talylor in 0 fino all’ordine 5 (a sinistra)
e fino all’ordine 7 (a destra), vediamo come l’approssimazione di sin(x) con
Pn (x) sia sempre migliore
al crescere di n.
§ 2.4
Esercizi risolti
1. Determinare la soluzione dell’equazione differenziale
y 0 = x2 y
tale che y(0) = 1.
Soluzione: Separando le variabili si trova
dy
= x2 dx ,
y
x3
lny =
+ c,
3
Z
1
dy =
y
y = ke
Z
x3
3
x2 dx ,
.
Imponendo y(0) = 1, si trova k = 1 e la funzione richiesta è y = e
§ 2.5
Altri esercizi
2. Data l’equazione differenziale
4y = xy 0 ,
S. Console – M. Roggero
x3
3
.
Equazioni differenziali
dire se
(a) y = x4 è una soluzione;
(b) trovare la soluzione generale.
3. Verificare che le funzioni y = x2 + cx sono soluzionidell’equazione differenziale
xy 0 − x2 − y = 0
e trovare la soluzione tale che y(1) = 0.
4. Determinare la soluzione dell’equazione differenziale
y0
=y
x
tale che y(0) = 1.
5. Determinare la soluzione dell’equazione differenziale
y0 =
x+1
.
y
6. Determinare la soluzione dell’equazione differenziale
y 0 = (x − 1)(y + 1)
tale che y(1) = 2.
Modelli Matematici applicati all’Ecologia (3.12.06)
15
Capitolo 3
Modelli di popolazioni isolate
Si definisce popolazione un insieme di esseri viventi, che condividono una porzione
di spazio, aventi delle caratteristiche comuni.
Generalmente si usa il termine popolazione per definire un gruppo di individui della
stessa specie, ma ha senso anche utilizzare questo termine individuando genericamente una popolazione di predatori, o per specificare all’interno della stessa specie
una popolazione di individui giovani, malati, maschi etc...
Una popolazione si puó considerare isolata se non è possibile individuare, all’interno dell’ecosistema, un partner privilegiato, una popolazione cioè la cui presenza
o assenza determini direttamente una mutazione sulla composizione della prima
popolazione.
Consideriamo popolazioni per cui è possibile fare la seguente ipotesi di omogeneità: in ogni istante ogni individuo può riprodursi o morire indipendentemente
dall’età.
In tal caso è possibile descrivere la variazione della popolazione in dipendenza
dal tempo e dalla densità della popolazione.
All’interno di simili popolazioni bisogna tenere conto della competetizione.
Infatti le forme e i modi con cui fertilità e mortalità possono essere alterate in
conseguenza dell’aumentata densità della popolazione possono essere molto diversi
a seconda delle popolazioni in esame. Un modo tradizionale di classificare i diversi meccanismi di competizione è quello di distinguere fra competizione per
interferenza e competizione per sfruttamento di risorse comuni.
La competizione per sfruttamento di risorse comuni è una competizione indiretta: gli individui possono al limite non incontrarsi mai, ma nonostante questo
si influenzano in maniera negativa. Il motivo di tale interazione risiede nel fatto
che, attraverso la propria semplice presenza, ogni individuo sottrae ad ogni altro
una parte delle risorse che sono essenziali alla sopravvivenza di entrambi. Animali
filtratori che competano per il particolato nella colonna d’acqua, oppure piante ad
alto fusto che con le loro chiome si ombreggino reciprocamente sottraendosi la luce
necessaria alla fotosintesi, sono solo due degli innumerevoli esempi possibili in tal
senso.
Il caso più semplice è quello di una singola popolazione isolata e in tal caso la
competizione è relativa all’uso delle risorse disponibili. L’effetto sarà più marcato al
16
Modelli di popolazioni isolate
17
crescere della densità di popolazione e, soprattutto, quando la densità della popolazione si avvicina alla soglia massima di sfruttamento delle risorse. Cominceremo
proprio con l’esaminare modelli di crescita di una popolazione isolata.
La competizione per interferenza è una competizione diretta: gli individui della
popolazione hanno contatti non mediati tra di loro, contatti che, nel caso animale,
si risolvono spesso in veri e propri scontri cruenti. Tale competizione dipende naturalmente dalla densità di individui presenti. Infatti, la probabilià che due individui si scontrino, e quindi si danneggino reciprocamente, cresce con l’abbondanza della popolazione. Nel capitolo 8 esamineremo modelli di crescita per popolazioni conviventi, tra cui il più antico e famoso è il modello preda-predatore di
Lotka-Volterra.
D’altra parte, un altro fattore di cui si deve tener conto è anche la cooperazione,
nel senso di interazione tra più popolazioni, in cui ciascuna beneficia della presenza
dell’altra. L’esempio più tipico in tal senso è la simbiosi. Prenderemo dunque in
esame anche modelli di cooperazione.
§ 3.1
Modelli di popolazioni isolate omogenee
Sia x(t) il numero di individui di una data popolazione in funzione del tempo t.
Sia ∆x la variazione di popolazione nel tempo ∆t.
L’ipotesi di omogeneità implica che la funzione che dà il numero dei nati nell’intervallo ∆t sia data da
n(x, t, ∆t) = ν(x, t)x(t)∆t ,
e, analogamente la funzione che dà il numero dei morti nell’intervallo ∆t sia data da
m(x, t, ∆t) = µ(x, t)x(t)∆t ,
Alternativamente, possiamo vedere ν(x, t)x(t)∆t e µ(x, t)x(t)∆t come primi termini
nello sviluppo di Taylor delle funzioni n(x, t, ∆t) e m(x, t, ∆t) in t = 0.
Le funzioni ν(x, t) e µ(x, t) sono rispettivamente il tasso istantaneo di natalità
e mortalità.
La differenza
r(x, t) = ν(x, t) − µ(xt, )
si dice tasso istantaneo di crescita.
Se non ci sono fenomeni migratori (popolazione isolata)
∆x = n(x, t, ∆t) − m(x, t, ∆t) = r(x, t)x(t)∆t ,
quindi
∆x
= r(x, t)x(t)
∆t
Modelli Matematici applicati all’Ecologia (3.12.06)
18
Capitolo 3
Se ora supponiamo ∆t infinitesimo, troviamo l’equazione differenziale che dà il
modello continuo:
dx
= r(x, t)x(t)
dt
Se invece supponiamo il tempo discreto, troviamo la successione ricorrente che
dà il modello discreto:
xn+1 − xn = r(xn )xn
x0
Vediamo ora modelli specifici che si possono ricavare per alcune scelte particolari
della funzione tasso istantaneo di crescita r(x, t). In generale, per costruire un
modello, possiamo procedere in due modi
ˆ scegliere la funzione r(x, t) che meglio si adatta alla particolari condizioni di
crescita;
ˆ scegliere a priori r(x, t) e studiare le corrispondenti curve di crescita.
Modello Malthusiano
La scelta più semplice è di supporre che il tasso di crescita sia costante. Questa
è quanto ha fatto Malthus alla fine del 1700 e per questo chiamiamo il modello
corrispondente modello Mathusiano.
Il modello continuo è dato dall’equazione differenziale
dx
= rx(t)
dt
(dove r è il tasso di crescita costante) che ha come soluzioni le curve esponenziali
x(t) = x0 ert
con x0 la popolazione iniziale.
Il modello discreto è la progressione geometrica
xn+1 = (1 + r)xn
x0
o, più esplicitamente
xn = x0 (1 + r)n
Notiamo che per r > 0 la popolazione raddoppia al tempo t̄ con
2x0 = x0 ert̄
e quindi per
t̄ =
log 2
.
r
S. Console – M. Roggero
Modelli di popolazioni isolate
19
Il tempo t̄ è detto tempo di raddoppio e posso scrivere
x(t) = x0 2t/t̄ .
Analogamente, per r < 0 la popolazione si dimezza per
t̄ = −
log 2
.
r
Il tempo t̄ è detto tempo di dimezzamento e posso scrivere
x(t) = x0 2−t/t̄ .
Il modello logistico
Il modello Malthusiano implica (nel caso di un tasso di crescita positivo) una crescita
esponenziale della popolazione. Per tale motivo non è spesso realistico, in quanto, se
la popolazione supera un certo livello, verranno a mancare le risorse e questo frenerà,
se non arresterà la crescita. In altre parole, si sentirà l’effetto della competizione
per sfruttamento di risorse comuni. Nel 1837 Verhulst propose una modifica al
modello Malthusiano, tenendo conto della carenza progressiva delle risorse al crescere
della popolazione.
Supponiamo che esista un valore massimo E detto popolazione di equilibrio
o capacità portante, superato il quale le risorse non siano più sufficienti e la
popolazione x(t) tenda a diminuire. Supponiamo poi che il tasso di crescita r(x)
sia proporzionale alla quantità delle risorse ancora disponibili, che è, a sua volta,
proporzionale a 1 − x/E. In altre parole
r(x) = k(1 − x/E) .
Troviamo allora l’equazione differenziale
x0 = k(1 − x/E)x = kx − k/Ex2
detta equazione di Verhulst o logistica.
Tale equazione si può risolvere separando le variabili. Infatti troviamo
dx
= kdt ,
x(1 − x/E)
Z
x
x0
Z
x
x0
dx
= kt ,
x(1 − x/E)
dx
+
x
Z
x
x0
1/Edx
= kt ,
1 − x/E
Modelli Matematici applicati all’Ecologia (3.12.06)
20
Capitolo 3
[log(x) − log(1 − x/E)]xx0 = kt ,
log
x(1 − x0 /E)
= kt .
x0 (1 − x/E)
Risolvendo questa equazione rispetto a x si trova la soluzione generale che è detta
funzione logistica
x(t) =
1.5
E
.
1 + (E/x0 − 1)e−kt
1.5
1
1
y(t)
y(t)
0.5
0
0.5
1
2
t
3
4
0
–0.5
–0.5
–1
–1
1
2
t
3
4
In figura viene rappresentato il campo di direzioni ed a fianco varie soluzioni al
variare della popolazione iniziale per l’equazione di Verhulst x0 = x(1 − x).
Osserviamo che
E
=E.
t→+∞ 1 + (E/x0 − 1)e−kt
lim
Dunque la popolazione di equilibrio E è l’ordinata dell’asintoto orizzontale destro
x = E.
Vedremo che x = E è un punto di equilibrio stabile del sistema.
S. Console – M. Roggero
Modelli di popolazioni isolate
21
Notiamo che dal confronto dei grafici della crescita malthusiana e logistica, si vede che il modello logistico è approssimato da quello malthusiano
quando la popolazione è lontana dalla capacità portante E. Se x << E il tasso di crescita per una
popolazione secondo il modello logistico è circa il
medesimo di quello malthusiano k. Diremo pertanto k il tasso di crescita “malthusiano”. Non appena
popolazione supera un certo livello, verranno a mancare le risorse e questo frenerà, se non arresterà la
crescita.
Un modo alternativo per arrivare ad un modello di crescita di tipo logistico è
il seguente. Consideriamo una popolazione animale con riproduzione continua la
cui dinamica, in condizioni di non affollamento, possa essere ben descritta da un
modello malthusiano continuo con tasso intrinseco di crescita positivo, pari ad r.
Supponiamo inoltre che:
ˆ gli animali si muovono in maniera casuale in un habitat di dimensione fissa;
ˆ la probabilità che due animali si incontrino è proporzionale al quadrato della
densità di animali presenti e quindi a x2 .
ˆ la probabilità che tre o più individui si incontrino contemporaneamente è
trascurabile (in altre parole trascuriamo i termini di ordine superiore nello
sviluppo di Taylor);
ˆ ogni volta che due individui si incontrano, tale incontro produce una diminuzione della natalità e/o un aumento della mortalità. In altre parole, un
aumento nella densità di individui presenti si traduce in una diminuzione del
tasso di crescita della popolazione.
Traducendo in equazioni le ipotesi appena avanzate, si ottiene l’equazione differenziale
x0 = rx − bx2
che si vede essere del tipo logistico.
Il coefficiente b è costante e positivo, detto coefficiente di competizione intraspecifica, che tiene conto in maniera aggregata sia della mobilità di ciascun individuo che
del danno prodotto da ciascun incontro. Confrontando l’equazione scritta ora con
l’equazione logistica, vediamo che la relazione tra il coefficiente di competizione
intraspecifica b e la popolazione di equilibrio E è
r
b= .
E
Modelli Matematici applicati all’Ecologia (3.12.06)
22
Capitolo 3
Possiamo anche considerare il modello discreto di crescita logistica.
Ricordando che
∆x
= r(x, t)x(t)
∆t
e che in questo caso si sceglie
r(x) = k(1 − x/E) ,
si trova la successione ricorrente
xn+1 = (1 + k(1 − xn /E))xn
x0
1
1
0.8
0.8
0.6
0.6
y
0.4
0.4
0.2
0.2
0
0.2
0.4
t
0.6
0.8
1
0
0.2
0.4
t
0.6
0.8
1
Scegliendo E = k = 1, visualizziamo i termini della successione disegnando la funzione y = x(2 − x) che genera per ricorrenza la successione e la bisettrice del 1
e 3 quadrante. Usando la diagonale si possono allora visualizzare i successivi termini (basta proiettare le linee verticali sull’asse x). Il diagramma ottenuto dà un
metodo estremamente semplice ed efficace per studiare la stabilità dei punti di equilibrio, introdotto negli anni ’50 dal statistico neozelandese Moran. Per tale ragione
chiameremo i diagrammi cosı̀ ottenuti diagrammi di Moran.
Si può visualizzare il comportamento del modello logistico per varie scelte della
popolazione di equilibrio e della popolazione iniziale, analizzando i diagrammi di
Moran delle successioni
xn+1 = a(1 − xn )xn
x0
al variare del parametro a con 0 ≤ a ≤ 4
Per 0 ≤ a ≤ 1 la successione tende alla popolazione di equilibrio nulla (estinzione).
S. Console – M. Roggero
Modelli di popolazioni isolate
23
Per 1 ≤ a ≤ 2 la successione tende alla popolazione di equilibrio (non nulla).
1
1
0.8
0.8
0.6
0.6
y
y
0.4
0.4
0.2
0.2
0
0.2
0.4
t
0.6
0.8
0
1
(a) a = 0.5, x0 = 0.8
0.2
0.4
t
0.6
0.8
1
(b) a = 1.5, x0 = 0.8
Per 2 < a ≤ 3 la successione tende sempre alla popolazione di equilibrio, ma
oscillando, con oscillazioni sempre più smorzate. Per 3 < a ≤ 4 ho invece un
comportamento caotico. Vedremo le ragioni di queste diverse situazioni nel capitolo
sulla stabilità.
1
1
0.8
0.8
0.6
0.6
y
y
0.4
0.4
0.2
0.2
0
0.2
0.4
t
0.6
0.8
(c) a = 2.7, x0 = 0.8
1
0
0.2
0.4
t
0.6
0.8
(d) a = 3.7, x0 = 0.8
Modelli Matematici applicati all’Ecologia (3.12.06)
1
24
Capitolo 3
§ 3.2
Modelli di popolazioni divise in due classi d’età
Consideriamo ora un modello in cui gli individui neonati non diventano fertili se non
dopo un certo periodo di tempo. Questo è più realistico per molte specie, incluso
ovviamente l’uomo.
Facciamo le seguenti ipotesi semplificatrici:
ˆ contiamo solo gli individui femmina, divisi in due gruppi d’età dette coorti di
ugual durata;
ˆ tutte le femmine di una stessa coorte, hanno la stessa fertilità e mortalità;
ˆ la popolazione è isolata;
ˆ le femmine della prima coorte non sono fertili e sono figlie di quelle della
seconda coorte;
ˆ superata l’età massima compresa nella seconda coorte, le femmine non sono
più fertili e quindi, anche se non muoiono, non saranno più considerate.
Indichiamo con xn e yn la successione temporale delle due classi d’età, dove
- xn denota la coorte non fertile;
- yn denota la coorte fertile.
Denotiamo inoltre con
`: la frazione della prima coorte che sopravvive ed entra a far parte della coorte
successiva delle adulte, ossia yn = `xn−1 ;
m: il numero medio di figli di ogni femmina fertile negli anni di riproduzione.
Le equazioni tempo discrete del modello sono:

 xn = myn

yn = `xn−1
1. caso: supponiamo m = m0 e ` = `0 costanti. Sostituendo la seconda equazione
nella prima, il modello si riduce ad un’unica successione ricorrente:
xn = m0 `0 xn−1 .
Il modello è dunque di tipo Malthusiano, essendo la successione geoemtrica di ragione
m0 `0 .
S. Console – M. Roggero
Modelli di popolazioni isolate
25
2. caso: supponiamo ora che le femmine nate in coorti numerose, facciano meno
figli di quelle nate in coorti esigue. Supponiamo sempre che ` = `0 sia costante, ma
ora
m = m0 (1 − yn /H) .
Le equazioni sono dunque ora

 xn = m0 (1 − yn /H)yn

yn = `0 xn−1
Sostituendo la seconda equazione nella prima, si elimina yn e si trova l’unica equazione
xn = m0 `0 (1 − `0 xn−1 /H)xn−1 .
Posto R0 = m0 `0 , si trova
xn = R0 xn−1 −
R0 `0 2
x
.
H n−1
Visualizziamo ora i termini della successione mediante diagrammi di Moran. La
funzione generatrice
R0 `0 2
y = R0 x −
x
H
incontra la bisettrice nell’origine e nel punto
xE =
H(R0 − 1)
.
R0 `0
Scegliendo R0 = 1.5 si trova, rispettivamente per termini iniziali x0 = 0.05 e x0 =
0.8, che i termini della successione convergono a xE .
1
1
0.8
0.8
0.6
0.6
y
y
0.4
0.4
0.2
0.2
0
0.2
0.4
t
0.6
0.8
1
0
0.2
0.4
t
0.6
Modelli Matematici applicati all’Ecologia (3.12.06)
0.8
1
26
Capitolo 3
Scegliendo invece R0 = 4.5 (e
x0 = 0.8 si ha invece un diagramma di Moran, che mostra una situazione di non convergenza ad xE .
2
1.8
1.6
1.4
1.2
y 1
0.8
0.6
0.4
0.2
0
§ 3.3
0.2
0.4
t
0.6
0.8
1
Le matrici
Per generalizzare quanto visto nel paragrafo precedente è necessario introdurre un
nuovo formalismo e metodo matematico fornito dalle matrici.
Chiamiamo matrici delle tabelle finite di elementi di un insieme N (in genere,
ma non sempre, un insieme di numeri) posti su ‘righe e colonne’. Nel seguito, fino a
diverso avviso N = R sarà sottinteso.
1 3 −2
è una matrice 2 righe e 3 colonne, brevemente 2 × 3
Esempio 3.1.
0 π
5
ad elementi in R.
Mm,n è l’insieme di tutte le matrici m × n ad elementi in N . Ogni matrice
A ∈ Mm,n si può scrivere come A = (aij ) con 1 ≤ i ≤ m, 1 ≤ j ≤ n e aij è l’elemto
che si trova all’incrocio della riga i e della colonna j.


a11 a12 a13 a14
Esempio 3.2. Per m = 3 e n = 4 si ha A =  a21 a22 a23 a24 .
a31 a32 a33 a34
Due matrici A = (aij ) ∈ Mmn e B = (bij ) ∈ Mm0 n0 coincidono se hanno le stesse
dimensioni m = m0 e n = n0 e gli elementi di posto corrispondente sono uguali ossia
se aij = bij per ogni coppia i, j.
1
Esempio 3.3. A = 1 2 e B =
sono diverse.
2
1 1
1 1
C=
eD=
sono diverse.
1 0
0 1
S. Console – M. Roggero
Modelli di popolazioni isolate
27
Le matrici di Mnn si dicono matrici quadrate. Se A = (aij ) è una matrice
quadrata, gli elementi aii che si trovano su righe e colonne con indice uguale formano
la diagonale principale di A.
Operazioni tra matrici
Somma di matrici La somma è definita solo tra matrici aventi le stesse dimensioni
e si esegue ‘posto per posto’.
In simboli, se A = (aij ) , B = (bij ) ∈ Mmn , la loro somma è data da:
A+B =C
dove
cij = (aij + bij )
La somma di matrici gode delle seguenti proprietà:
1. Associativa e commutativa.
2. Esiste l’elemento neutro rispetto alla somma ossia la matrice nulla 0 di
Mm,n (R) che ha 0 in ogni posto ed è tale che A + 0 = 0.
3. Esiste l’opposto di ogni matrice A = (aij ) ossia una matrice −A tale che
A + (−A) = 0 che è −A = (−aij ).
Nota bene: il simbolo 0 può indicare tante matrici di dimensioni diverse, tutte
fatte da soli zeri. Non usiamo simboli diversi per distinguerle perchè di volta in volta
sarà chiaro dal contesto quale stiamo considerando.
Prodotto per uno scalare Siano A una matrice e λ un numero reale. Col simbolo
λA indicheremo la matrice B, con le stesse dimensioni di A, ottenuta moltiplicando
per λ ogni elemento di A:
λA = λ(aij ) = (λaij )
Diremo che due matrici non nulle con le stesse dimensioni A e B sono proporzionali se B = λA per un qualche λ ∈ R.
Prodotto righe per colonne Consideriamo innanzi tutto una ‘matrice riga’
 e una

b1
 b2 

‘matrice colonna’ della stessa lunghezza:
R = (a1 a2 . . . an ) e C = 
 .. .
bn
Il loro prodotto RC è il numero a1 b1 + a2 b2 · · · + an bn .
Se A = (aij ) ∈ Mm,n e B = (bjk ) ∈ Mm0 ,n0 sono due matrici di dimensioni tali
che n = m0 , definiamo il loro prodotto A · B come la matrice P di dimensioni m × n0
in cui l’elemento di posto ik è il prodotto della riga i-esima di A e della colonna
k-esima di B ossia
pik = ai1 b1k + ai2 b2k + · · · + aim bmk .
Modelli Matematici applicati all’Ecologia (3.12.06)
28
Capitolo 3
Esempio 3.4.
=
1 2 3
−1 0 5


π 0 2
·  −3 1 1  =
1 1 0
1 · π + 2 · (−3) + 3 · 1
1·0+2·1+3·1
1·2+2·1+3·0
−1 · π + 0 · (−3) + 5 · 1 −1 · 0 + 0 · 1 + 5 · 1 −1 · 2 + 0 · 1 + 5 · 0
=
π−3 5
4
−π + 5 5 −12
=
Nota bene: La somma è definita solo tra matrici con le stese dimensioni, mentre
il prodotto è definito solo se le righe della prima matrice sono lunghe come le colonne
della seconda.
Il prodotto righe per colonne non è commutativo. Infatti:
ˆ può darsi che A · B sia definito, ma B · A no;
ˆ se anche A·B e B ·A sono entrambe definite, possono avere dimensioni diverse;
ˆ se anche A · B e B · A sono entrambi definite e hanno le stesse dimensioni
(questo accade se A, B sono matrici quadrate con le stesse dimensioni), i due
prodotti possono comunque essere matrici diverse.
1 4
0 2
1 1
mentre B · A =
, A·B =
B=
Esempio 3.5. A =
1 6
1 2
2 1
4 2
.
5 3
Nel seguito ometteremo il simbolo · per il prodotto.
Anche se non è commutativo, valgono comunque le seguenti proprietà del prodotto:
1. Proprietà associativa: se A, B, C sono matrici tali che AB e BC sono
definite, allora si può fare il prodotto delle tre e vale l’uguaglianza A (BC) =
(AB) C.
2. Proprietà distributive: se A, B, C sono matrici tali che A (B + C) è definito, allora si ha A (B + C) = AB +AC (e analogamente (B + C) A = BA+CA
ogni volta che le operazioni sono definite).
3. Esistenza dell’identità moltiplicativa: se In è la matrice n × n che ha
elementi tutti nulli, tranne quelli sulla diagonale principale che valgono 1,
allora ogni volta che il prodotto è possibile si ha AIn = A e In B = B.
S. Console – M. Roggero
Modelli di popolazioni isolate
Esempio 3.6. Le matrici identià 2 × 2 e 3 × 3 sono

1
1 0
I2 =
,
I3 =  0
0 1
0
29
rispettivamente :

0 0
1 0 .
0 1
Nel seguito scriveremo semplicemente I omettendo l’indece n quando la dimesnione sarà chiara dal contesto.
1 1
0 −1
Esempio 3.7. Consideriamo le matrici A =
e B =
. Se
−1 0
1
1
eseguiamo i prodotti nei due modi possibili troviamo AB = BA = I.
Se A e B sono matrici quadrate n×n e soddisfano le condizioni AB = I e BA = I,
diciamo che B è l’inversa di A (e A è l’inversa di B): in tal caso si scrive B = A−1 .
Importanti differenze rispetto alle operazioni tra numeri sono mostrate dagli
esempi seguenti.
1 1
1
1
Esempio 3.8. Consideriamo le due matrici A =
e B =
,
1 1
−1 −1
che sono
diverse dalla matrice nulla 0. Eseguendo il prodotto si ottiene:
entrambe
0 0
= 0.
AB =
0 0
Diciamo quindi che tra le matrici non vale la legge di annullamento del
prodotto.
Esempio
le matrici A e B dell’esempio precedente e la matrice
3.9. Consideriamo
3 −2
C =
, che sono tutte diverse da 0. Eseguendo i prodotti si ottiene:
−3
2
AB = AC (perchè in entrambi i casi si trova 0), anche se B 6= C.
Diciamo quindi che tra le matrici non vale la legge di cancellazione.
0 0
(che è diversa da 0) e
Esempio 3.10. Consideriamo la matrice D =
0 1
0 0
a b
un’altra matrice qualsiasi X =
. Il prodotto delle due è: DX =
c d
c d
che non può in ogni caso essere la matrice identità.
Quindi non tutte le matrici hanno l’inversa. In generale soltanto le matrici
quadrate possono avere l’inversa ed in tal caso ne hanno una sola.
a b
Esempio 3.11. La matrice 2 × 2 generica del tipo X =
ha inversa se e
c d
soltanto se il numero ∆ = ad − bc non è nullo.
Si può infatti verificare con calcoli
d/∆ −b/∆
diretti che se ∆ 6= 0, allora la matrice Y =
è l’inversa di X.
−c/∆
a/∆
Modelli Matematici applicati all’Ecologia (3.12.06)
30
Capitolo 3
Il determinante di una matrice quadrata
Il numero ∆ che abbiamo introdotto alla fine del paragrafo precedente si chiama
determinante della matrice X. Più in generale, ad ogni matrice quadrata A ∈ Mn,n
si può associare un numero che indicheremo con det(A), il suo determinante, con
la proprietà:
la matrice A ha inversa se e soltanto se det(A) 6= 0.
Il determinante può essere definito in modo intrinseco; noi ci limiteremo però a
dire come può essere calcolato. Vedremo alcuni modi diversi per calcolare det(A): è
importante sottolineare il fatto, tutt’altro che evidente, che tutti questi modi diversi
portano allo stesso risultato.
Sviluppo del determinante rispetto a una riga Sia A = (aij ). Fissiamo una
riga qualsiasi di A: sia la riga i = k. Allora:
det(A) = (−1)k+1 ak1 det(Ak1 ) − ak2 det(Ak2 ) + · · · + (−1)n−1 akn det(Akn )
dove la somma è a segni alterni e Akj è la matrice (n − 1) × (n − 1) che si ottiene
da A ‘cancellando’ la riga k-esima e la colonna j-esima.
Si ripete quindi il procedimento con le sottomatrici ottenute, fino ad arrivare a
sottomatrici 2 × 2 il cui determinante si calcola direttamente con la formula vista.
Esempio 3.12. Sviluppiamo il seguente determinante rispetto alla terza riga:
0
1
det @ 3
1
2
0
1
1
»
„
1
2
−1 A = (−1)3+1 1 · det
0
2
1
−1
«
„
− 1 · det
1
3
1
−1
«
„
+ 2 · det
1
3
2
0
«–
=
= (−2 − 0) − (−1 − 3) + 2(0 − 6) = −10
Sviluppo del determinante rispetto a una colonna Sia A = (aij ). Fissiamo
una colonna qualsiasi di A: sia la colonna j = h. Allora:
det(A) = (−1)1+h a1h det(A1h ) − a2h det(A2h ) + · · · + (−1)n−1 anh det(Anh )
dove la somma è a segni alterni e Aih è, come prima, la matrice (n − 1) × (n − 1)
che si ottiene da A ‘cancellando’ la riga i-esima e la colonna h-esima.
Si ripete quindi il procedimento con le sottomatrici ottenute, fino ad arrivare a
sottomatrici 2 × 2.
Esempio 3.13. Sviluppiamo il determinante della matrice A dell’esempio precedente rispetto alla seconda colonna:
3 −1
1 1
1
1
1+2
det(A) = (−1)
2 · det
− 0 · det
+ 1 · det
=
1
2
1 2
3 −1
= −2(6 + 1) + 0 − 1(−1 − 3) = −10.
S. Console – M. Roggero
Modelli di popolazioni isolate
31
Come si può osservare dal confronto di questi due metodi, il determinante non
cambia se si scambiano ordinatamente le righe di una matrice con le sue colonne,
ottenendo la cosiddetta matrice trasposta t A:




1
3 1
1 2
1
0 1 .
Esempio 3.14. A =  3 0 −1  t A =  2
1 −1 2
1 1
2
Proprietà del determinante Il determinante di una matrice gode di molte importanti proprietà. Ne elenchiamo alcune. Su di esse si basa il metodo migliore per
calcolare i determinanti, conveniente soprattutto quando le dimensioni sono grandi,
che sarà spiegato nel prossimo paragrafo.
1. Se una matrice ha una riga oppure una colonna tutta di zeri, il suo
determinante è 0. (Basta sviluppare rispetto alla riga o alla colonna tutta
nulla).
2. Se si moltiplicano per una stessa costante λ tutti gli elementi di una
riga o di una colonna, il determinante risulta moltiplicato per λ.
3. Se si scambiano tra loro due righe (oppure due colonne) il determinante cambia segno.
4. Se una matrice ha due righe oppure due colonne uguali (o proporzionali), allora il determinante è 0. Infatti scambiando tra loro le
righe (o le colonne) uguali, il determinante deve cambiare segno (per il punto
3) ) ma anche rimanere invariato (perchè la matrice rimane invariata). L’unico
numero uguale al suo opposto è 0.
Osserviamo che nel calcolo di un determinante conviene sempre scegliere lo
sviluppo rispetto ad una riga o ad una colonna che contengano tanti zeri, cosı̀ da
ridurre al minimo il numero di determinanti di sottomatrici che si devono calcolare.
Ancora più utile è ridurre la matrice cioè trasformarla in un’altra che abbia lo
stesso determinante, ma che contenga molti più elementi nulli. Le proprietà precedenti permettono appunto di operare trasformazioni su una matrice senza alterare
il suo determinante o almeno con alterazioni conosciute a priori. Questo metodo va
sotto il nome di riduzione della matrice.
Riduzione di una matrice.
Questo metodo può essere adoperato sempre, anche nel caso di matrici non quadrate,
poichè ha molteplici applicazioni oltre a quella del calcolo dei deteminanti.
Consideriamo una generica matrice A anche non quadrata. Operiamo su di essa
mediante le seguenti trasformazioni sulle righe:
Modelli Matematici applicati all’Ecologia (3.12.06)
32
Capitolo 3
i) moltiplicare tutti gli elementi di una riga per una stessa costante λ;
ii) scambiare tra loro due righe;
iii) sommare ad ogni elemento di una riga il corrispondente elemento di un’altra
fissata riga moltiplicato per una costante λ.
Una opportuna sequenza di queste trasformazione sulle righe di una matrice A ∈
Mm,n permette di ottenere una matrice A0 ∈ Mm,n ridotta ossia una matrice A0
tale che:
in ogni riga che non sia tutta nulla vi è un elemento non nullo al di sotto del
quale vi sono solo zeri.
Si parla di matrice ridotta a scalini se la matrice è ridotta e gli elementi
speciali di ciascuna riga si incontrano procedendo da sinistra verso destra e dall’alto
verso il basso.
Esempio 3.15. Esempio di matrice

1
 0

 0
0
ridotta con elementi speciali evidenziati:

2 0
3
5 0 −1 
.
2 0
0 
0 0
0

1
 0
Esempio di matrice ridotta a scalini: 
 0
0
2
0
0
0

0
3
2 −4 
.
0 −3 
0
0
Algoritmo per trasformare una matrice in una a scalini
ˆ Si cerca la prima colonna non tutta nulla.
ˆ Si scambiano se necessario le righe in modo che in quella colonna l’elemento
più in alto sia non nullo.
ˆ Si aggiunge ad ogni riga successiva un opportuno multiplo della prima in modo
da ottenere zeri al di sotto di quell’elemento.
ˆ Si procede allo stesso modo sulle righe e colonne successive, senza più alterare
o usare le righe sovrastanti.
Quando il procedimento termina si ottiene una matrice a scalini.
Applicazioni:
1. Il determinante di una matrice quadrata a scalini è il prodotto degli elementi
della diagonale principale.
S. Console – M. Roggero
Modelli di popolazioni isolate
33
2. Per calcolare il determinante di una matrice si procede alla sua riduzione,
tenendo conto che il determinante viene moltiplicato per λ ogni volta che si
moltiplica per λ una riga e cambia segno ogni volta che due righe vengono
scambiate tra loro.
3. Per calcolare l’inversa di una matrice A ∈ Mn,n (con determinante non nullo) si
affianca alla matrice A la matrice In e si opera sulle righe con le trasformazioni
i), ii) e iii) in modo che al posto di A si ottenga In . La matrice n × n che si
trova dove prima c’era In è l’inversa di A.
4. Per risolvere i sistemi lineari o i sistemi di equazioni differenziali lineari, si
riduce la matrice dei coefficienti.
Sistemi lineari
Un sistema lineare di m equazioni in n incognite:

a11 x1 + a12 x2 + · · · + a1n xn = b1



a21 x1 + a22 x2 + · · · + a2n xn = b2
......



am1 x1 + am2 x2 + · · · + amn xn = bm
può essere pensato in modo matriciale come AX = B dove A ∈ Mm,n è la matrice
dei coefficienti, X ∈ Mn,1 è la colonna delle incognite e B ∈ Mm,1 è la colonna dei
termini noti:






a11
a12 . . . a1n
x1
b1
 a21




a22 . . . a2n 
 X =  x2  B =  b2 
A=
 ......

 ... 
 ... 
am1 am2 . . . amn
xn
bm
Notiamo che se A è una matrice invertibile, il sistema AX = B avrà una sola
soluzione data da:
X = A−1 B.
Se in particolare il sistema lineare è omogeneo, ossia la matrice dei termini noti
B è tutta nulla, allora AX = 0 avrà la sola soluzione nulla X = 0 quando A ha
inversa e quindi quando det(A) 6= 0, mentre avrà anche soluzioni non nulle quando
det(A) = 0.
Sia ora A|B la matrice ottenuta affiancando ad A una ulteriore colonna data da
B. Le operazioni di riduzione su A|B non alterano le soluzioni del sistema.
Possiamo allora ricondurci al caso AX = B con A|B ridotta.
Il sistema ha soluzioni se e soltanto se le matrici A e A|B hanno lo stesso numero
r di righe non tutte nulle. In tal caso le soluzioni dipendono da n − r parametri
liberi.
Modelli Matematici applicati all’Ecologia (3.12.06)
34
Capitolo 3
Autovalori e autovettori di una matrice quadrata
Da ora in poi chiameremo vettore colonna o brevemente vettore una matrice con
una sola colonna.
Sia A una matrice quadrata n × n e sia v un vettore n × 1. Diciamo che v
è un autovettore della matrice A se il prodotto Av è un vettore proporzionale a
v. In formule: Av = λv. Se esiste un tale vettore non nullo, diremo che λ è un
autovalore di A.
In altri termini un autovettore è una soluzione del sistema lineare AX = λX.
Portando tutte le indeterminate a primo membro, otteniamo un sistema lineare
omogeneo con matrice dei coefficienti A0 = A − λIn .
Condizione necessaria e sufficiente perchè λ sia autovalore di A è allora:
det(A − λIn ) 6= 0.
Le soluzioni dell’equazione det(A − tIn ) = 0 nell’incognita t sono quindi gli
autovalori di A; se λ è una di tali soluzioni, gli autovettori relativi si determinano
risolvendo il sistema lineare omogeneo (A − λIn )X = 0.
1 −5
. Allora
Esempio 3.16. Sia A =
−3
3
1 − t −5
= t2 − 4t − 12.
det(A − tI2 ) = det
−3 3 − t
Le soluzioni di det(A − tI2 ) = 0 sono λ1 = −2 e λ2 = 6. Gli autovettorirelativi
5k
all’autovalore −2, ossia le soluzioni X di (A + 2I2 )X = 0, sono i vettori
3k
con k ∈RR; gli
autovettori
relativi
a
6
ossia
le
soluzioni
di
(A
−
6I
)X
=
0
sono
i
2
h
vettori
con h ∈ R.
−h
Talvolta il polinomio det(A − tIn ) = 0 non ha soluzioni in R. Conviene allora
calcolare le sue soluzioni complesse ossia nell’insieme dei numeri complessi C.
Potenze di una matrice
Nelle applicazioni può essere necessario calcolare potenze An di una matrice quadrata
A con un esponente n molto alto. Non conviene eseguire gli n prodotti successivi
AA, A2 A, . . . , An−1 A, ma il calcolo può essere notevolmente più breve se si ricorre
a uno dei seguenti risultati.
Teorema 3.17. (Cayley-Hamilton) Se F (t) è il polinomio caratteristico della
matrice quadrata A, ossia il polinomio det(A − tIn ), allora F (A) = 0.
S. Console – M. Roggero
Modelli di popolazioni isolate
35
Per calcolare An si divide il polinomio tn per F (t) ottenendo tn = F (t)Q(t)+R(t)
dove il polinomio resto R(t) ha grado ≤ n − 1. Quindi An = R(A).
Teorema 3.18. Se A ∈ Mn,n possiede n autovalori diversi (reali oppure complessi),
allora esiste una matrice P tale che A = P DP −1 , dove D è una matrice diagonale
ossia ha tutti elementi nulli tranne quelli sulla diagonale principale.
La matrice D ha sulla diagonale principale gli n autovalori distinti di A; la
matrice P ha come colonne un autovettore per ognuno degli n autovalori di A nello
stesso ordine con cui compaiono in D.
Si ha quindi An = P Dn P −1 dove Dn si ottiene semplicemente elevando alla
potenza n-esima gli elementi della diagonale di D.
§ 3.4 Modelli di popolazioni strutturate per età o per
taglia
Possiamo a questo punto generalizzare quanto visto nel paragrafo sui modelli di
popolazioni divise in due classi d’età, considerando un numero qualsiasi r di coorti
ottenute dividendo in r periodi di pari durata il ciclo vitale (medio) degli individui
di una certa popolazione.
Conviene considerare come unità di misura per il tempo la durata di un periodo.
La dinamica della popolazione strutturata in fasce d’età nel tempo è data dalla
sequenza dei vettori vn di lunghezza r aventi come i-esima entrata il numero xi (n)
di femmine che all’istante t = n si trova a far parte della i-esima coorte.
Indichiamo con ai la percentuale di femmine nella coorte i che ragiunge la coorte
successiva (per le posizioni fatte all’inizio del §3.2 si ha ar = 0) e con bi il numero
medio di figlie nate durante un periodo da ogni femmina della i − esima fascia d’età.
Sia infine A la matrice r × r che ha come prima riga (b1 , . . . , br e i-esima riga
(successiva alla prima) tutta nulla tranne ai nella posizione i − 1. Vale allora la
relazione
Avn = vn+1 .
La matrice A viene chiamata matrice di Leslie.
Supponendo che i dati della matrice A siano costanti nel tempo, gli autovettori
con autovalore 1 (purché a componenti tutte positive!) sono le distribuzioni per
fasce d’età che permettono alla popolazione di rimanere invariata.
Gli autovettori con autovalore λ > 0 (purché a componenti tutte positive!) sono
le distribuzioni per fasce d’età che, pur modificandosi numericamente, mantengono
inalterate le distribuzioni percentuali per fascia d’età. Se λ > 1 la popolazione
cresce, se 0 < λ < 1, la popolazione decresce.
Modelli Matematici applicati all’Ecologia (3.12.06)
36
Capitolo 3
Esempio 3.19. Una certa popolazione animale ha un ciclo vitale che in condizioni
ottimali dura circa 15 anni.
Suddividiamo la popolazioni (considerando solo le femmine) in tre coorti per
fasce d’età: da 0 a 5, da 5 a 10 e da 10 a 15. Indichiamo con x1 (t), x2 (t), x3 (t) la
consistenza numerica delle tre fasce d’età nell’istante t. Sia t = 0 l’istante iniziale
di osservazione e t = 1 indichi un periodo di 5 anni.
Indichiamo poi con ai la percentuale di femmine nella coorte i = 1, 2, 3 che
raggiungono la coorte successiva (ovviamente a3 = 0) e con bi il numero medio di
figlie nate durante un periodo da ogni femmina nella i − esima fascia d’età.
Le consistenze numeriche dei tre gruppi nel momento t = n sono legate a quelle
nel momento t = n + 1 dalle relazioni:

 x1 (n + 1) = b1 x1 (n) + b2 x2 (n) + b3 x3 (n)
x2 (n + 1) = a1 x1 (n)

x3 (n + 1) = a2 x2 (n)
Possiamo esprimere sinteticamente questa relazione mediante la notazione matriciale
X(n + 1) = AX(n) dove X(t) è il vettore colonna con le tre componenti xi (t) e


b1 b 2 b 3
A =  a1 0 0  .
0 a2 0
§ 3.5
Esercizi
1. Siano date le matrici:
1
A=
0
2
3
B=
−1
1
2 −2

0 −1 1
−1
2 1  E=
D= 2
0
1
0 1

1 −1
2
2
C=
1
0
1
0
0
3


−1 −1
F =  0 −3 
−2
1
i) Eseguire tutte le possibili operazioni tra le matrici A, . . . , F .
ii) Calcolare il determinante delle matrici A, B e D sviluppando rispetto a righe e colonne
e mediante la riduzione.
iii) Dire quali delle precedenti matrici hanno inversa e calcolare tali inverse quandoi
esistono.
iv) Calcolare autovalori e autovettori reali delle matrici A, B e D.
S. Console – M. Roggero
Modelli di popolazioni isolate
2. Ridurre la matrice:
37


0
2
0 0 0
 1 −1
1 0 2 
.
G=
 3
1 −2 2 1 
4
0 −1 2 3
3. Risolvere i seguenti sistemi lineari:

 2x − 3y + z = −1
x − 3y − z = −2

2x + y + +3z = −2

 2x + 3y + z − t = 0
y − z + 2t = 1

2x + y + 3z − 5t = −2

x+y =0




 y+z =0
z+t=0


x+y+z+t=0



x+z =0
4. Una certa popolazione umana è suddivisa in rurale e urbana: indichiamo con Rt e Ut
la consistenza numerica dei due gruppi all’istante t (in milioni di individui). Supponiamo
di esaminare la popolazioni in istanti successivi discreti, ad esempio sia t = 0 il momento
iniziale dell’osservazione e sia poi t = 1 un anno.
Le consistenze numeriche dei due gruppi nell’istante t = n sono legate a quelle nell’istante
t = n + 1 dalle relazioni:
Rn+1 = (a − b(1 − c))Rn + bcUn
Un+1 = b(1 − c)Rn + (a − bc)Un
Il significato delle costanti a, b, c è il seguente:
ˆ a = fattore di crescita di entrambe le popolazioni;
ˆ c = rapporto ottimale tra la popolazione rurale e quella totale ai fini dell’approvvigionamento alimentare;
ˆ b = tasso di migrazione dalla campagna verso la città relativo all’eccesso di popolazione
rurale rispetto a quella ottimale.
Supposto che i valori delle costanti non varino nel tempo e che siano a = 1.2, b = 1.1 e
c = 1/3 studiare l’evoluzione delle due popolazioni nei tre casi seguenti:
R0
20
R0
2
R0
5
;
=
.
=
;
=
U0
10
U0
30
U0
10
5. Considerando la situazione illustrata nell’Esempio 3.19 Calcolare l’evoluzione di una
popolazione con x1 = x2 = x3 = 1000 con matrice


0
4 3
0 0 .
A =  1/2
0 1/4 0
Cosa succede dopo 50 anni? e dopo 500?
Vi sono distribuzioni che si mantengano costanti nel tempo?
Modelli Matematici applicati all’Ecologia (3.12.06)
38
Capitolo 3
Vi sono distribuzioni che si mantengano percentualmente costanti nel tempo?
Rispondere alle stesse domande nel caso la matrice A sia:


0
0 6
0 0 .
A =  1/2
0 1/3 0
§ 3.6
Spunti per approfondimenti
Altri modelli per popolazioni isolate
Uno dei modelli più utilizzati per descrivere la dinamica di popolazioni a riproduzione concentrata con dipendenza da densità è stato introdotto da Beverton e
Holt (1957) con riferimento agli stock ittici. L’idea sulla quale si sono basati è
data dal seguente schema tratto dalle note di Marino Gatto e Renato Casagrandi
reperibili alla pagina web http://olmo.elet.polimi.it/ecologia/dispensa/
Le notazioni sono le seguenti:
ˆ La quantità Nt rappresenta il numero di adulti (o di femmine adulte se la
popolazione è sessuata) che costituiscono la popolazione della stagione t-esima.
ˆ f indica la fertilità di ogni adulto.
ˆ Il numero di uova (o di uova da cui nasceranno femmine) complessivamente
presenti nella popolazione è pari a Et = f Nt .
Dopo la deposizione delle uova tutti gli adulti riproduttori muoiono come nel caso
della cavalletta. Non tutte le uova deposte dagli individui si schiudono con successo
e non tutti gli avannotti da esse originati sopravvivono fino allo stadio di giovani
pesci dando origine a giovani individui.
In figura è rappresentata la sopravvivenza da uovo a giovane pesce σt in funzione
del numero di uova deposte E
S. Console – M. Roggero
Modelli di popolazioni isolate
39
La sopravvivenza delle uova nel modello di
Beverton-Holt è supposta decrescere con la densità delle uova stesse.
Il parametro σ0 rappresenta la sopravvivenza dallo stadio di uovo a quello di
giovane in condizioni di scarso affollamento, mentre il parametro ρ regola l’entità
della dipendenza da densità. In seguito, solo una frazione σY degli Yt giovani presenti
riuscirà ad emergere come adulto alla stagione successiva. Si trova allora il modello
tempo-discreto
σY σ0 f N t
Nt+1 = σY Yt = σY σE Et =
.
1 + ρf Nt
Un altro modello è quello logistico con ritardo: è basato sulle seguenti discrepanze tra modello logistico e realtà naturale:
ˆ Secondo il modello logistico l’effetto della densità sul tasso di accrescimento
è istantaneo. Nella realtà c’è invece un ritardo tra variazione di densità ed
effetto sul tasso di crescita r.
ˆ Le condizioni ambientali e quindi la capacità portante E sono supposti costanti
per il modello logistico. Nella realtà spesso E può variare con la stagione o
negli anni.
ˆ Si suppone che tutti gli individui abbiano le stesse caratteristiche, mentre in
una popolazione si riscontrano individui con caratteristiche differenti.
Un modello logisto con ritardo è descritto dall’equazione differenziale
x0 = k(1 − x(t − τ )/E)x(t) ,
dove τ è il ritardo di tempo.
Un esempio è dato dal seguente grafico per la crescita di una popolazione di
pecore.
Modelli Matematici applicati all’Ecologia (3.12.06)
40
Capitolo 3
Accenniamo infine brevemente al modello logistico con prelievo costante.
A rigore non si tratta più di una popolazione isolata, ma l’interazione con l’ambiente
è solo il prelievo (ad esempio quello da parte dell’uomo in una riserva).
È descritto dall’equazione differenziale
x0 = k(1 − x(t)/E)x(t) − H ,
dove H è il prelievo costante (“harvesting” in inglese). Tale equazione è di difficile soluzione ma è possibile analizzarla qualitativamente, ad esempio considerando
l’analogo modello tempo discreto e costruendo il diagramma di Moran. Vedremo
in un prossimo capitolo il suo esame in relazione alla stabilità delle soluzioni di
equilibrio.
Tabelle e grafi di vita e matrici di Leslie
Per modelli di popolazioni divise in classi di età sono importanti le tabelle di vita
o life tables, di cui è riportato un esempio.
S. Console – M. Roggero
Modelli di popolazioni isolate
41
Le colonne riportano informazioni relative alla mortalità e fertilità specifica per
ogni classe o età.
Il ciclo di vita di una popolazione divisa in classi di
età è schematizzato dal cosiddetto grafo di vita.
Vediamolo su un esempio di popolazione a riproduzione preriodica. Supponiamo
che il ciclo sia di questo tipo:
ˆ la vita è di al più di 3 anni;
ˆ le femmine sono in ugual numero dei maschi;
ˆ le femmine del primo anno non sono fertili, quelle di due anni 8 piccoli, le
femmine di tre anni, 6 piccoli.
ˆ il 40% dei piccoli sopravvive dalla nascita al primo anno, l’80% dal primo
nascita al secondo anno e il 70% dal secondo al terzo.
Indichiamo con fi il tasso finito di fertilità degli individui della coorte i-esima.
Nell’esempio:
f1 = 0 ,
f2 = 8 ,
f3 = 6 .
Denotiamo con si la frazione di femmine della coorte i-esima che sopravvivono sino
all’età i + 1. Nel caso dell’esempio:
s0 = 0.4 ,
s1 = 0.8 ,
s2 = 0.7 .
Il grafo di vita è composto da nodi e archi:
ˆ Il nodo i-esimo ( ni ) rappresenta il numero di femmine che prima dell’inizio
della stagione riproduttiva di un determinato anno t hanno età i.
ˆ Gli archi a tratto continuo indicano il contributo delle madri in termini di nuove
nate, mentre gli archi tratteggiati rappresentano il processo di invecchiamento.
Modelli Matematici applicati all’Ecologia (3.12.06)
42
Capitolo 3
L’esempio precedente, per una
popolazione a riproduzione periodica divisa in 3 classi di età dà
luogo al grafo di vita riportato in
figura.
Le informazioni contenute nelle tabelle e grafo di vita permettono di costruire
la matrice di Leslie. Indichiamo con ni (t) la popolazione (femminile) della coorte iesima al tempo t. Nell’esempio avremo allora tre funzioni n1 (t), n2 (t), n3 (t) (ovvero
un vettore n(t) = (n1 (t), n2 (t), n3 (t))) cosı̀ determinate:
Invecchiamento L’equazione che dà le femmine che hanno più di un anno di età
(i > 1) è codificata negli archi tratteggiati:
ni+1 (t + 1) = si ni (t) .
Nell’esempio
n2 (t + 1) = 0.8 n1 (t) ,
n3 (t + 1) = 0.7 n2 (t) .
Riproduzione L’equazione che regola la riproduzione è determinata dagli archi
continui e quindi dai termini si fj . Nell’esempio
n1 (t + 1) = 0.4 · 4n2 (t) + 0.4 · 3n2 (t) = 1.6n2 (t) + 1.2n2 (t) .
Le consistenze numeriche dei tre
momento t + 1 dalle relazioni:

 n1 (t + 1)
n2 (t + 1)

n3 (t + 1)
gruppi nel momento t sono legate a quelle nel
= 1.6n2 (t) + 1.2n3 (t) ,
= 0.8 n1 (t) ,
= 0.7 n2 (t) .
Dunque la matrice di Leslie è


0 1.6 1.2
0
0 .
A =  0.8
0 0.7
0
S. Console – M. Roggero
Modelli di popolazioni isolate
43
In generale



A=


s0 f1 s0 f2
s1
0
0
s2
...
...
0
0
...
...
...
...
...
s0 fm−1 s0 fm
0
0
0
0
...
...
sm−1
0






Per approfondimenti su entrambi gli argomenti cui si è accennato in questo paragrafo, si consultino testi [BC, F], e si vedano anche le slides del seminario “Modello
logistico con ritardo e popolazioni strutturate per età” alla pagina web del corso.
Modelli Matematici applicati all’Ecologia (3.12.06)
Capitolo 4
Calibrazione di modelli matematici
Supponiamo che siano disponibili conteggi o stime di una data popolazione in
stagioni successive.
Ad esempio, consideriamo i dati per la quantità di piante della specie phlox
drummondii in una località del Texas (dati raccolti da Leverich e Levine nel 1979)
tempo (mesi)
piante
1
996
2
158
3
154
4
151
5
147
6
136
7
105
8
74
9
22
I dati dell’andamento della popolazione sono nella colonna “Nx ”.
Il nostro scopo è:
ˆ decidere quale modello matematico sia il più adatto a descrivere l’andamento
della popolazione;
ˆ una volta deciso il modello, stimare nel modo più preciso i parametri del modello. Questa stima dei parametri viene detta calibrazione di dati o “fitting”
del modello.
Vediamo come procere nel caso dei dati sulla phlox drummondii.
Iniziamo con il riportare i dati su un foglio
Excel. Non si tratta di un programma professionale, ma tuttavia è di uso abbastanza agevole
oltre che estrema diffusione.
44
Calibrazione di modelli matematici
45
Possiamo usare il comando ‘‘grafico’’
(‘‘chart’’ nella versione in inglese) per
fare un grafico della distribuzione di popolazione. Dal grafico si vede come la distribuzione sembri approssimativamente rispettare il modello malthusiano, cioè la descrescita sia di tipo esponenziale.
Ricordiamo che la una legge di tipo malthusiano ha equazione
y = aebt .
Per vedere se l’andamento è di tipo esponenziale “linearizziamo”: calcoliamo il
logaritmo di ambo i membri dell’equazione precedente trovando
log y = log a + bt .
Perciò in un riferimento semilogaritmico (ascissa t e ordinata log(y)) l’equazione
diventa lineare.
L’idea è allora di calcolare (nella colonna C della tabella di Excel) log(y). A tale
scopo scriviamo della cella C3 la formula
=LN(B3).
Selezionando le altre colonne e scegliendo il comando ‘‘riempimento in basso’’
(‘‘fill down’’ nella versione in inglese)
dal menu ‘‘Modifica’’ (‘‘Edit’’) si possono inserire automaticamente i dati per le
altre colonne (cioè automaticamente viene
inserita la formula
=LN(B4)
nella quarta riga etc.)
Modelli Matematici applicati all’Ecologia (3.12.06)
46
Capitolo 4
Disegnamo ora il grafico con ‘‘chart’’
inserendo come asse delle ascisse il tempo e asse delle ordinate i dati della colonna C. Al tempo stesso inseriamo la retta
di regressione ‘‘trendline’’ ( ‘‘linea
di tendenza’’ nella versione italiana) insieme con la sua equazione.
Ricordiamo dal corso di statistica che,
dato un insieme di punti nel piano, la retta
di regressione lineare rappresenta la retta
che “passa più vicino ai punti”. In termini più precisi, viene minimizzata la somma
dei quadrati delle distanze dei punti stessi
dalla retta; il grafico fornisce un’intuizione
del procedimento.
Nel caso della distribuzione di phlox drummondii la retta di regressione lineare
ha equazione
y = −0, 3066x + 6, 4189
Il coefficiente di correlazione, che in questo caso è
R2 = 0, 7261
fornisce una misura di quanto i dati seguono un andamento lineare: 0 ≤ R2 ≤ 1, il
valore 1 è assunto se i dati sono perfettamente allineati e tanto più R2 si avvicina a
1, tanto più i dati si avvicinano a seguire un andamento lineare.
Possiamo concludere che nel nostro caso il modello malthusiano non è ottimale
(R2 non è molto vicino a 1), ma se vogliamo usare tale modello la scelta migliore
per i parametri a e b nell’equazione
y = aebt
S. Console – M. Roggero
Calibrazione di modelli matematici
47
è data da
log a = 6, 4189 ,
b = −0, 3066 .
Calcoliamo ora a come EXP(6,4189) e determiniamo i dati secondo il modello
malthusiano: scriviamo una nuova colonna con la formula:
= $H$3*EXP($F$4*(A3))
che ricopiamo in tutta la colonna D usando il comando ‘‘fill down’’. Si inserisce
il nome della cella tra dollari (ad esempio $H$3) affinchè se uno ricopia automaticamente i dati delle celle tale dato non sia aumentato ogni volta di uno come da
‘‘fill down’’.
Disegnamo ora il grafico dei dati sperimentali in confronto con quelli del modello
malthusiano calibrato. Si noti che c’è una
certa discrepanza, come del resto prevedibile, dato che il coefficiente di correlazione
non è molto vicino a 1.
Modelli Matematici applicati all’Ecologia (3.12.06)
48
Capitolo 4
§ 4.1
Trasformazioni linearizzanti
Il procedimento usato per il modello malthusiano può essere applicato in altri casi.
Il primo passo è di linearizzare il modello usando trasformazioni linearizzanti
Tabella 1: Alcuni modelli e trasformazioni linearizzanti
Modello
Esponenziale
Logistico
Sigmoide
Equazione
y = aebt
y = K/(1 + ae−bt )
y = K/(1 + atb )
Linearizzazione
log y = log a + bt
log(K/y − 1) = log a − bt
log(K/y − 1) = log a + b log t
La linearizzazione per il modello esponenziale (malthusiano) si ottiene facilmente
prendendo il logaritmo di ambo i membri dell’equazione y = aebt . In altre parole, ottengo un grafico lineare in un riferimento semilogaritmico (ascisse t, ordinate
log(y)).
Per il modello logistico, la linearizzazione si ottiene con i seguenti passaggi
algebrici:
y =
1
y
1
y
K
y
=
=
K
,
1 + ae−bt
1 + ae−bt
,
K
1 + ae−bt
,
K
= 1 + ae−bt ,
K
− 1 = ae−bt ,
y
K
log
−1
= log(ae−bt ) ,
y
K
log
−1
= log a − bt .
y
I calcoli per il modello sismoide sono molto simili.
Vediamo ora come si può calibrare un modello di crescita logistica.
§ 4.2
Calibrazione per un modello di crescita logistica
Partiamo dai seguenti dati per una popolazione di parameci aurelia che inseriamo
in un foglio elettronico di Excel
S. Console – M. Roggero
Calibrazione di modelli matematici
49
Possiamo usare il comando ‘‘chart’’
di Excel per fare un grafico della distribuzione di popolazione. Dal grafico si vede
come la distribuzione sembri rispettare il
modello logistico.
Come prima operazione, vogliamo effettuare un fitting con il modello logistico di
crescita mediante trasformazione delle variabili (in modo a ridursi a una equazione
lineare) e poi mediante regressione lineare.
Ricordiamo che l’equazione logistica
y = K/(1 + ae−bt )
diventa lineare con la trasformazione linearizzante
log(K/y − 1) = log a − bt .
Dai dati supponiamo che la popolazione limite (di equilibrio) sia K = 560 (ricordiamoci che è l’asintoto orizzantale della funzione logistica). Determiniamo allora a
e b mediante regressione lineare.
L’idea è allora di calcolare (nella colonna C della tabella di Excel) log(K/y − 1).
A tale scopo scriviamo della cella C2 la formula =LN(560/B2-1). Otteniamo:
Modelli Matematici applicati all’Ecologia (3.12.06)
50
Capitolo 4
Disegnamo ora il grafico con ‘‘chart’’
inserendo come asse x i giorni e asse y i
dati della colonna C. Al tempo stesso inseriamo la retta di regressione ‘‘trendline’’
insieme con la sua equazione.
Dunque
log(a) = 4, 332 ,
a = 76, 0963 ,
b = 0, 6143 .
Si trova la curva logistica
y = 560 ∗ 1/(1 + 76.0963 ∗ exp(−.6134 ∗ t)).
In alternativa si puo’ effettuare mediante lo strumento ‘‘analisi di dati’’
(nella versione in inglese si trova sotto il menu ‘‘Tools’’ e si chiama ‘‘Data
analysis’’.)
Osservazione 4.1. Potrebbe succedere che non si trovi tale opzione tra gli strumenti. Bisogna in primo luogo cercare tra gli ‘‘Add-ins’’ (‘‘Aggiunte’’) e -se lo
si trova- aggiungere “analisi di dati”. Potrebbe pero’ non trovarsi. In questo caso
tale aggiunta non e’ stata installata (con l’installazione standard di ‘‘Office’’ non
viene installata): pertanto e’ necessario usare il CD di installazione.
Vediamo come possiamo ottimizzare la scelta dei parametri usando il ‘‘Solver’’
(‘‘Risolutore’’) di Excel. (Come nel caso di ‘‘analisi di dati’’, potrebbe
S. Console – M. Roggero
Calibrazione di modelli matematici
51
succedere che non si trovi tale opzione tra gli strumenti. Si proceda come nella nota
precedente.)
Nella seconda colonna inseriamo i dati in base alla formula data dal modello di
crescita logistica y = K/(1 + aexp(−bt)). Per fare questo, si inserisce ad esempio
nella terza colonna la formula scrivendo ad esempio nella 2 riga (relativa a 2 giorni)
= $H$2/(1+($H$3*EXP(-$H$4*A2))). Notare che che uso il parametro K inserito
nella cella H2, a inserito nella cella H3, etc. Come le altre volte, se si inserisce il
nome della cella tra dollari (come in $H$2), usando (‘‘fill down’’) il dato non
viene cambiato.
Successivamente calcoliamo l’errore quadratico inserendo in ogni riga della colonna D formule come =(C2-B2)ˆ2 (questo per la cella D2) etc. Otteniamo:
Usando il “Solver” (“risolutore”) di Excel (che si trova sotto il menu ‘‘Tools’’
Modelli Matematici applicati all’Ecologia (3.12.06)
52
Capitolo 4
- ‘‘Strumenti’’) ottimizziamo i parametri K, a e b minimizzando la somma degli
scarti quadratici medi.
Bisogna impostare il “Solver” indicando che i parametri sono nelle celle e H2-H4
(e quindi inserendo il comando $H$2:$H$4 sotto ‘‘changing cells’’) degli scarti quadratici medie la somma degli scarti quadratici medi D16 (‘‘Target cell’’
D16). Si trova:
Possiamo infine fare il grafico dei dati misurati e di quelli calibrati.
S. Console – M. Roggero
Calibrazione di modelli matematici
§ 4.3
Esercizi
1. Ripetere la procedura fatta per i parameci aurelia usando il modello sigmoide.
2. Determinare la retta di regressione per la crescita del bacillo aspergillus niger
3. Analizzare in modo analogo
Modelli Matematici applicati all’Ecologia (3.12.06)
53
54
Capitolo 4
la crescita del saccharomyces cerevisiae come nella tabella
la crescita di girasoli come nella tabella
S. Console – M. Roggero
Capitolo 5
Sistemi di equazioni differenziali
Molti problemi sono governati non da una singola equazione differenziale, ma da
un sistema di più equazioni. Ad esempio questo succede se si vuole descrivere un
sistema ecologico di due o più popolazioni.
Consideriamo quindi un sistema di due equazioni differenziali
x0 = f (x, y)
y 0 = g(x, y)
(2)
dove x = x(t) e y = y(t) sono due funzioni incognite (ad esempio le densità di due
popolazioni al variare del tempo) e f (x, y) e g(x, y) sono funzioni assegnate (che
supponiamo differenziabili in una regione del piano). Un sistema di questo tipo si
dice autonomo poichè le funzioni f (x, y) e g(x, y) non dipendono dal tempo t.
Esempio 5.1.
x0 = −y
y0 = x
è un sistema detto oscillatore armonico. Più in generale,
x0 = −5x + 2y
y 0 = 2x − 2y
è un sistema lineare omogeneo a coefficienti costanti. È utile rappresentare
il sistema in forma matriciale:
x0
y0
=
−5 2
2 −2
x
y
Una soluzione del sistema (2) è data da una coppia di funzioni x(t), y(t) derivabili
in un certo intervallo. In generale, le soluzioni del sistema non sono uniche ma
dipendono da due costanti arbitrarie. Per determinarle, purchè valgano condizioni
di regolarità per le funzioni f e g, è sufficiente fissare le condizioni iniziali, Vale cioè
il
55
56
Capitolo 5
Teorema 5.2 (di Cauchy). Sotto opportune ipotesi di regolarità per le funzioni f
e g, assegnate le condizioni iniziali
x(t0 ) = x0 ,
y(t0 ) = y0 ,
esiste un’unica coppia di funzioni x(t), y(t), soluzioni del sistema (2) che verificano
le condizioni iniziali.
Esempio 5.3 (L’oscillatore armonico). Abbiamo visto che un esempio di sistema
omogeneo a coefficienti costanti è l’oscillatore armonico:
0
x = −y
y0 = x .
Ora, se deriviamo la prima equazione, troviamo
x00 = −y 0
e sostituendo la seconda equazione
x00 = −x .
Quest’ultima è un’equazione del secondo ordine a coefficienti costanti che ha come
soluzione generale
x = c1 cos t + c2 sin t .
Sostituendo nella seconda equazione si trova
y 0 = c1 cos t + c2 sin t
e quindi y = c1 sin t − c2 cos t. Non è difficile vedere che le curve soluzione sono
circonferenze. Ad esempio se la condizione iniziale è il passaggio per il punto P0 =
(1, 0) si ha
c1 = 1 ,
c2 = −1
e quindi la soluzione
x = cos t
y = sin t .
che è la circonferenza di centro l’origine e raggio unitario.
Si noti che al variare del parametro t, ogni soluzione x = x(t), y = y(t) descrive
nel piano x, y una curva che viene detta orbita oppure traiettoria del sistema.
Un punto (x0 , y0 ) tale che f (x0 , y0 ) = 0 e g(x0 , y0 ) = 0 si dice punto stazionario
o di equilibrio. Si noti che la soluzione dell’equazione differenziale passante per
(x0 , y0 ) è (x(t), y(t)) = (x0 , y0 ). In altre parole, l’orbita di un punto stazionario è il
punto stesso.
S. Console – M. Roggero
Sistemi di equazioni differenziali
§ 5.1
57
Il piano delle fasi
Un sistema di equazioni differenziali
x0 = f (x, y)
y 0 = g(x, y)
si può interpretare in modo “fisico” nel modo seguente:
la coppia di funzioni (f (x, y), g(x, y)) definisce in ogni punto del piano un vettore
e quindi descrive un campo vettoriale V . Le soluzioni (x(t), y(t)) possono essere
interpretate come posizione di una particella al tempo t la cui velocità (x0 , y 0 ) è data
in ogni punto dal campo V . Le orbite sono quindi le traiettorie della particella. Il
piano x, y viene detto piano delle fasi.
Esempio 5.4. Nel caso del sistema
x0 = x + 2y
y 0 = 3x + 2y
che può essere scritto in forma matriciale come
0 x
1 2
x
=
0
y
3 2
y
scriviamo alcuni vettori del campo:
−1
=⇒
P =
1
P =
P =
−3
−2
2
0
V =
1 2
3 2
−1
1
=
1
−1
=⇒
=⇒
1 2
2
2
V =
=
3 2
0
6
−7
1 2
−3
=
V =
−13
3 2
−2
Se facciamo questo per un gran numero di vettori otteniamo il piano delle fasi.
Il seguente disegno rappresenta il piano delle fasi e alcune traiettorie.
Modelli Matematici applicati all’Ecologia (3.12.06)
58
Capitolo 5
Esempio 5.5. Nel caso dell’oscillatore armonico
x0 = −y
y0 = x
abbiamo il seguente piano delle fasi.
S. Console – M. Roggero
Sistemi di equazioni differenziali
59
Le orbite (eccetto l’orbita per l’origine, che è l’origine stessa) sono circonferenze
di centro l’origine.
Se consideriamo il sistema (che possiamo pensare come una lieve pertubazione
dell’oscillatore armonico)
0
x = −0.05x − y
y 0 = x − 0.05y
abbiamo il seguente piano delle fasi con orbite che a differenza del caso dell’oscillatore armonico non sono chiuse.
Modelli Matematici applicati all’Ecologia (3.12.06)
60
Capitolo 5
§ 5.2 Soluzioni di sistemi lineari omogenei a coefficienti
costanti
Vediamo ora come si possono trovare in generale le soluzioni di un sistema lineare
omogeneo a coefficienti costanti
0
x = ax + by
y 0 = cx + dy
che può essere scritto in forma matriciale come
0 a b
x
x
=
c d
y
y0
ossia:
Z 0 = AZ ,
posto
Z=
x
y
,
0
Z =
x0
y0
,
A=
a b
c d
.
L’idea è di operare in modo simile al caso delle equazioni del second’ordine omogenee
a coefficienti costanti: cerchiamo soluzioni di “tipo esponenziale”, cioè del tipo
Z = W eλt
S. Console – M. Roggero
Sistemi di equazioni differenziali
61
con W un vettore che scriviamo come matrice 2 × 1 (vettore colonna)
w1
W =
,
w2
e λ = α + iβ un numero complesso.
Sostituendo Z = W eλt nell’equazione Z 0 = AZ troviamo
λW eλt = AW eλt
cioè, dividendo per eλt :
λW = AW .
Questo significa che λ è un autovalore di A e W un autovettore di A. Dunque
la ricerca delle soluzioni del sistema lineare omogeneo Z 0 = AZ si riduce alla ricerca
degli autovalori ed autovettori della matrice A.
Esempio 5.6. Consideriamo il sistema
0 x
−5 2
x
.
=
0
y
2 −2
y
Autovalori ed autovettori di A sono:
1
2
2
−1
λ1 = −1 ,
W1 =
λ2 = −6 ,
W2 =
.
Ogni coppia autovalore-autovettore produce una soluzione:
1
2
e−6t .
e−t ,
2
−1
La soluzione generale è una combinazione lineare delle soluzioni:
1
2
−t
e + c2
e−6t ,
W = c1
2
−1
con c1 e c2 costanti arbitrarie.
Leggendo il risultato “per righe” possiamo scrivere le soluzioni generali per x e y:
x = c1 e−t + 2c2 e−6t ,
y = 2c1 e−t − c2 e−6t .
In generale, abbiamo i seguenti casi per le soluzioni di sistemi lineari omogenei
a coefficienti costanti:
Modelli Matematici applicati all’Ecologia (3.12.06)
62
Capitolo 5
ˆ A ha autovalori reali distinti :
la soluzione generale è del tipo
W = c1 W1 eλ1 t + c2 W2 eλ2 t ,
dove λ1 e λ2 sono gli autovalori, W1 e W2 i relativi autovettori e c1 e c2 costanti
arbitrarie.
ˆ A ha autovalori reali coincidenti :
in tal caso possiamo avere due casi a seconda che in corrispondenza all’unico
autovalore λ di molteplicità 2, A possieda 2 autovettori indipendenti oppure
questo non avvenga. Il secondo caso è più complesso e non lo trattiamo rimandando ad esempio a [Antonio C. Capelo, Modelli matematici in biologia,
introduzione all’ecologia matematica, Zanichelli Decibel editore, Cap. I.5.1].
Se invece A ha 2 autovettori indipendenti W1 e W2 , la soluzione generale è del
tipo
W = c1 W1 eλt + c2 W2 eλt ,
dove c1 e c2 sono costanti arbitrarie.
ˆ A ha autovalori complessi coniugati λ1 = α + iβ, λ2 = α − iβ:
in questo caso gli autovettori sono W1 ± iW2 e usando la formula di Eulero si
vede che la soluzione generale è del tipo
W = c1 eαt (W1 cos(βt) − W2 sin(βt)) + c2 eαt (W1 sin(βt) + W2 cos(βt)) ,
dove c1 e c2 sono costanti arbitrarie.
Esempio 5.7. Consideriamo il sistema
0 0 1
x
x
.
=
0
−4 0
y
y
Gli autovalori della matrice del sistema A sono soluzioni dell’equazione caratteristica:
−λ 1 2
−4 −λ = λ + 4 = 0
e sono quindi dati da ±2i. I corrispondenti autovettori sono
1
±2i
cioè
1
0
±i
0
2
.
S. Console – M. Roggero
Sistemi di equazioni differenziali
63
La soluzione generale è quindi:
W
1
0
cos(2t) −
sin(2t) +
= c1
0 2
1
0
sin(2t) +
cos(2t) ,
+c2
0
2
con c1 e c2 costanti arbitrarie.
Leggendo il risultato “per righe” possiamo scrivere le soluzioni generali per x e y:
x = c1 cos(2t) + c2 sin(2t) ,
y = −2c1 sin(2t) + 2c2 cos(2t) .
Notiamo che possiamo scrivere le equazioni delle orbite in modo “più espressivo”
quadrando ed eliminando il tempo t tra le equazioni:
x2 = c21 cos2 (2t) + 2c1 c2 sin(2t) cos(2t) + c22 sin2 (2t)
y 2 = 4c21 sin2 (2t) − 8c1 c2 sin(2t) cos(2t) + 4c22 cos2 (2t)
Moltiplicando la prima equazione per 4 e sommando le due equazioni, troviamo
4x2 + y 2 = 4(c21 + c22 ) ,
che rappresenta l’equazione di un ellisse. Dunque, in questo caso le orbite sono
ellissi, come si vede dalla seguente figura, che mostra il piano delle fasi e alcune
orbite.
Modelli Matematici applicati all’Ecologia (3.12.06)
64
§ 5.3
Capitolo 5
Esercizi
S. Console – M. Roggero
Capitolo 6
Stabilità
§ 6.1
Punti di equilibrio di un’equazione differenziale
Consideriamo ora equazioni differenziali o sistemi di equazioni differenziali autonomi.
Una equazione differenziale in forma normale del primo ordine è autonoma se
è della forma
x0 = F (x)
cioè la funzione F è indipendente dal tempo t.
Un sistema di equazioni differenziali

dx



 dt = f (x, y)



 dy = g(x, y)
dt
è autonomo se le funzioni f e g sono indipendenti da t.
Cominciamo dal caso di una equazione differenziale. Un punto x0 tale che
F (x0 ) = 0 si dice punto stazionario o di equilibrio. Si noti che la soluzione
dell’equazione differenziale passante per x0 è x(t) = x0 .
1.5
1
y(t)
0.5
0
–0.5
–1
1
2
t
3
4
Esempio 6.1 (Modello logistico). Nel caso
dell’equazione logisitica x0 = x(1 − x), i punti
di equilibrio sono x = 0 e x = 1; si noti che
il punto 1 è la popolazione di equilibrio (nel
caso generale è il punto x = E). Dalla figura
si vede che le soluzioni passanti per i punti
con x = 0 e x = 1 sono le funzioni costanti,
rispettivamente x(t) = 0 e x(t) = 1.
65
66
Capitolo 6
Cosa succede se perturbo leggermente un sistema che si trova in stato di equilibrio? ci sono essenzialmente due tipi di equilibrio:
instabile e stabile.
Equilibrio: Stabile o instabile?
Nel caso del modello logistico, abbiamo visto che ci sono due punti di equilibrio: x =
0 e x0 = E. Vedremo che uno dei due punti è di equilibrio stabile, l’altro instabile.
Intuitivamente è un sistema stabile quello tipo il quadro in figura, instabile, tipo la
matita. In questa parte introduttiva non daremo ancora una definizione formale di
stabilità, che vedremo solo nel paragrafo sui sistemi di equazioni differenziali.
Un punto x0 è di equilibrio
stabile se il sistema ritorna
allo stato di equilibrio dopo
una piccola perturbazione. Un
punto x0 è di equilibrio instabile se il sistema si allontana dallo stato di equilibrio
dopo una piccola perturbazione. In alcuni testi si usa il termine asintoticamente stabile invece di stabile in quanto
l’equilibrio è raggiunto quando t → ∞.
Osservazione 6.2. Questi due casi di stabilità ed instabilità sono i più semplici.
Come detto, torneremo sul concetto di equilibrio paragrafo sui sistemi di equazioni
differenziali. Notiamo che nei casi esaminati i punti sono isolati. Se invece i punti
sono di accumulazione di punti di equilibrio possono presentarsi altre situazioni
più complesse. In generale, si parla di punti di equilibrio semistabile, metastabile,
stabile-metastabile, instabile-metastabile. Un punto è di equilibrio semistabile se il
sistema ritorna all’equilibrio o se ne allontana ulteriormente a seconda di come viene
perturbato. Se invece il sistema non ritorna allo stato di equilibrio originale, ma
non si allontana troppo da esso, assumendo immediatamente o asintoticamente una
configurazione di equilibrio vicina alla precedente si dice metastabile. In questi casi
i punti stazionari sono isolati. Invece i punti di equistabile-metastabile, instabilemetastabile non sono isolati, e non li prenderemo in considerazione - si veda ad
esempio [C, Cap. I.6.2] per la loro definizione e maggiori dettagli.
Stabilità per i punti di equilibrio di una equazione differenziale
Esaminiamo la stabilità dei 2 punti di equilibrio del modello logistico.
S. Console – M. Roggero
Stabilità
67
dx
> 0
dt
e quindi la popolazione cresce (il
punto nel grafico si muove verso
destra). Se x < 0 oppure x > E
(ovviamente, x < 0 non ha senso dal punto di vista biologico),
la popolazione decresce (il punto
sul grafico si muove verso sinistra). Le frecce mostrano che il
punto di equilibrio x = 0 è instabile, mentre il punto x = E è
stabile.
Se 0 < x < E, allora
Dal punto di vista naturalistico, questo significa che dopo una piccola variazione
dell’entità di popolazione da x = 0 (ad esempio, a causa dell’immigrazione di un
piccolo numero di individui), la popolazione non torna mai indietro a questo punto
di equilibrio. Invece, la popolazione cresce fino a quando raggiunge (all’infinito) lo
stato di equilibrio x = E. Una deviazione dallo stato x = E torna a questo stato di
equilibrio. In generale, sia data l’equazione differenziale
dx
= f (x) .
dt
Ciò che differenzia i punti di equilibrio stabile e instabile è il coefficiente angolare
della retta tangente alla funzione y = f (x) nel punto di equilibrio:
ˆ se è positivo, il punto è di equilibrio instabile;
ˆ se è negativo, il punto è di equilibrio stabile.
Osservazione 6.3. Potrebbe pero’ succedere che il coefficiente angolare nel punto
di equilibrio sia nullo (cioè dx/dt = 0). In questo caso è necessario effettuare lo
sviluppo di Taylor della funzione y = f (x) di ordine superiore al primo (quanto
fatto prima è in altre parole uno sviluppo di Taylor del primo ordine). Per maggiori
dettagli si veda ad esempio [C, Cap. I.6.2].
Stabilità e biforcazione
Le equazioni differenziali che descrivono un fenomeno come la crescita di una popolazione contengano uno o più parametri. La domanda che ci poniamo ora è la
seguente: come dipendono dai parametri gli equilibri?
Esempio 6.4 (modello logistico). Per il modello logistico
x0 = k(1 − x(t)/E)x(t)
Modelli Matematici applicati all’Ecologia (3.12.06)
68
Capitolo 6
abbiamo due parametri: E che è la capacità portante del sistema e k che rappresenta
il tasso di crescita “malthusiano” (cioè finchè il modello è lontano dalla capacità
portante). Si può ridurre il sistema ad un solo parametro misurando la popolazione
relativa alla capacità portante E: con il cambiamento di variabile u := Ex l’equazione
diventa
u0 = k(1 − u(t))u(t)
(equivalentemente si può pensare di di fissare E = 1.) Analizzeremo nell’Esempio 6.8 la dipendenza dal parametro nel caso dell’analogo modello tempo-discreto.
Esempio 6.5 (modello logistico con prelievo costante). Supponiamo ora di
avere un prelievo costante H (“harvesting”)
x0 = k(1 − x(t)/E)x(t) − H ,
dove H è il prelievo costante. Analogamente a prima, misuriamo la popolazione
relativa alla capacità portante (cioè poniamo u := Ex ) e inoltre normalizzaimo tra
loro H e k. Posto h := H
k , consideriamo cioè l’equazione ad un unico parametro:
u0 = (1 − u(t))u(t) + h .
I punti di equilibrio sono le
soluzioni dell’equazione
u2 − u + h = 0
e sono quindi dati da
√
1 ± 1 − 4h
∗
u =
2
Pertanto come si vede in figura
si sono due punti di equilibrio u∗1
e u∗2 per h < 14 , questi due punti
coincidono per h = 14 e cessano di esistere punti di equilibrio per h > 14 . Si osservi
che ovviamente se h = 0 ci riduciamo al modello logistico, che ha due punti equilibrio inquesto caso: 0 è instabile, mentre l’altro punto (che coincide con la capacità
portante) è di equilibrio stabile.
S. Console – M. Roggero
Stabilità
69
Questa situazione cambia di poco se
il prelievo è piccolo rispetto alla capacità
portante, cioè per h < 41 . I due punti di
equilibrio u∗1 < u∗2 sono rispettivamente di
equilibrio instabile e stabile.
Quindi nell’esempio sul modello logistico con prelievo costante, quando
il prelievo è piccolo, ci sono 2 equilibri di cui uno instabile e l’altro stabile.
Quando il prelievo aumenta oltre la soglia h = 14 l’equilibrio scompare.
Diciamo che c’è una biforcazione per
h = 41 .
Questo perchè h = 14 è il valore per cui
c’è un cambiamento significativo nel comportamento della popolazione.
Si può disegnare un diagramma detto diagramma di biforcazione in cui
si rappresentano le popolazioni di equilibrio u∗ al variare di h (si noti che
non è una funzione ad un valore!). Nel
caso dell’esempio il diagramma è una
parabola con l’asse parallelo a quello
delle ascisse. Il ramo superiore corrisponde agli equilibri stabili u∗2 .
Modelli Matematici applicati all’Ecologia (3.12.06)
70
Capitolo 6
§ 6.2
Stabilità per gli equilibri di modelli discreti
Consideriamo ora un sistema tempo-discreto
xn+1 = f (xn ) .
La successione ricorrente è dunque descritta dalla funzione y = f (x).
Il punto p è un punto fisso della funzione f se f(p) = p. I punti fissi di f sono i
punti fissi o di equilibrio del sistema tempo-discreto.
Ci chiediamo ora quando i punti fissi o di equilibrio sono stabili o instabili. In
questo caso si parla anche di
ˆ attrattori (punti di equilibrio stabile)
ˆ repulsori (punti di equilibrio instabile).
L’idea è sempre quella di prendere un punto q vicino ad uno fisso. Consideriamo
la successione ricorrente
xn+1 = f (xn ) .
con valore iniziale x0 = q. Se la successione converge a p (cioè i valori di xn per
n grande si avvicinano sempre più a p) il punto p è un attrattore. Altrimenti se i
valori di xn si allontano sempre più da p, il punto p è un repulsore.
Esempio 6.6 (Modello logistico). Nel caso della crescita logistica rappresentata
dalla successione ricorrente
xn+1 = (2 − xn ))xn
x0
i punti fissi sono p = 0 e p = 1 (p = E nel caso generale). Ci si può rendere conto
del diverso comportamento per i 2 punti fissi o di equilibrio usando una delle tante
Applet che si trovano in rete (ad esempio
http://www.lboro.ac.uk/departments/ma/gallery/doubling/ oppure anche dalle
seguenti figure:
S. Console – M. Roggero
Stabilità
71
caso del punto fisso p = 1: stabile o attrattore.
caso del punto fisso p=0 instabile o repulsore.
Come si possono riconoscere i punti di equilibrio stabile (attrattori) o instabile
(repulsori) a partire da proprietà della funzione f ?
Si può enunciare, analogamente al caso continuo (equazioni differenziali) un
semplice criterio in base al coefficiente angolare della setta tangente nel punto fisso.
Consideriamo un modello tempo-discreto: (ora la variabile è il numero naturale
n). Sia p un punto fisso.
Il punto p è un attrattore se
Modelli Matematici applicati all’Ecologia (3.12.06)
72
Capitolo 6
df −1 <
<1
dx p
df dove è
il coefficiente angolare della retta tangente in p alla curva y = f (x).
dx p
df Altrimenti, se
è maggiore di 1 o minore di −1, p è un repulsore.
dx p
Ad esempio, nel caso della crescita logistica (f (x) = x(2 − x)), nel punto 1 il
coefficiente angolare della retta tangente è 0 (< 1) e il punto 1 è effettivamente un
attrattore.
Esempio 6.7 (Un modello di popolazione divisa in due classi d’età). Abbiamo visto (pagina 24) che un modello tempo-discreto di popolazioni in due classi
d’età xn e yn è regolato, sotto determinate ipotesi dalle equazioni:

 xn = m0 (1 − yn /H)yn

yn = `0 xn−1
Eliminando yn e si trova l’unica equazione
xn = R0 xn−1 −
R0 `0 2
x
,
H n−1
S. Console – M. Roggero
Stabilità
73
dove si è posto R0 = m0 `0 . Analizziamo ora la stabilità del punto di equilibrio
xE =
H(R0 − 1)
.
R0 `0
La derivata della funzione generatrice f (x) = R0 x −
R 0 `0 2
H x
in xE è data da
f 0 (xE ) = −R0 + 2 .
Dunque si ha equilibrio stabile se
−1 < −R0 + 2 < 1
cioè
1 < R0 < 3 .
Questo è in accordo con le situazioni mostrate nelle figure a pagina 25.
Esempio 6.8 (modelli logistici al variare di un parametro). Consideriamo
ora la famiglia di modelli logistici tempo discreti già considerata a pagina 22.
xn+1 = fa (xn ) = a(1 − xn )xn
x0
al variare del parametro a con 0 ≤ a ≤ 4
Se x0 = 0 o x0 = 1, oppure a = 0, la successione è costante 0 (per n ≥ 1).
Altrimenti il comportamento della successione dipende fortemente dal valore di a.
Vedremo che il comportamento è regolare (indipendente dalla popolazione iniziale)
per 0 ≤ a ≤ 3 e caotico (o dipendente dal valore iniziale) per 3 < a ≤ 4.
ˆ Per 0 ≤ a ≤ 1, l’unico punto di equilibrio nel primo quadrante è 0 e la successione è monotona decrescente e tende alla popolazione di equilibrio nulla
(estinzione). D’altra parte fa0 (0) = a, quindi 0 ≤ fa0 (0) ≤ 1 e la stabilità segue
dal criterio.
ˆ Per 1 < a ≤ 2 la successione tende (in modo monotono crescente per n ≥ n0 )
alla popolazione di equilibrio
1
`=1− .
a
È facile vedere che in questo caso la stabilità segue dal criterio. Infatti:
fa0 (`) = 2 − a
e quindi in questo caso 0 ≤ fa0 (`) < 1.
Modelli Matematici applicati all’Ecologia (3.12.06)
74
Capitolo 6
1
1
0.8
0.8
0.6
0.6
y
y
0.4
0.4
0.2
0.2
0
0.2
0.4
t
0.6
0.8
0
1
(a) a = 0.5, x0 = 0.8
0.2
0.4
t
0.6
0.8
1
(b) a = 1.5, x0 = 0.8
1
1
0.8
0.8
0.6
0.6
y
y
0.4
0.4
0.2
0.2
0
0.2
0.4
t
0.6
0.8
(c) a = 2.7, x0 = 0.8
1
0
0.2
0.4
t
0.6
0.8
(d) a = 3.7, x0 = 0.8
S. Console – M. Roggero
1
Stabilità
75
ˆ Per 2 < a ≤ 3 la successione tende sempre alla popolazione di equilibrio `, ma
oscillando, con oscillazioni sempre più smorzate. Si noti che −1 ≤ fa0 (`) < 0
in questo caso. Per a = 3 non si può allora applicare il criterio, ma si vede
direttamente che ` è stabile.
I punti stazionari della successione corrispondono a orbite 1−periodiche. per
trovare le orbite 2−periodiche consideriamo i punti fissi della successione ottenuta
iterando fa , cioè
xn+1 = fa ◦ fa (xn ) = a2 xn (1 − xn )(ax2n − axn + 1) .
Ovviamente 0 e ` sono punti fissi di questa successione. Gli altri punti fissi si trovano
risolvendo l’equazione
1
ax2 − (1 + a)x + 1 + = 0 .
a
Il discriminante di questa equazione è ∆ = (1 + a)(a − 3), quindi per a > 3 abbiamo
altri 2 punti fissi.
Per a = 3, ∆ = 0 e l’equazione diventa 3(3x − 2)2 = 0 e la seconda iterata ha
sempre un’unica orbita stabile 23 . Per a > 3 quest’orbita diventa instabile e i due
punti fissi addizionali producono orbite 2-periodiche stabili. Il comportamento della
successione diventa sempre più complicato all’avvicinarsi di a a 4.
Disegnamo ora l’appicazione che a
ogni a ∈ [0, 4] associa le orbite periodiche stabili della successione logistica. Si può vedere che per a ∈ [0, 3]
questa è un’applicazione liscia
a tratti:
è infatti data da a 7→ max 0, 1 − a1 .
Tuttavia appena a supera il valore 3,
l’applicazione diventa a più valori e il
suo grafo effettua una serie di biforcazioni sempre più complicate.
Per ulteriori dettagli rimandiamo all’articolo di Albert Milani “Dynamical systems: regularity and chaos”, Rendiconti del Seminario Matematico dell’Università e del Politecnico di Torino (Vol 64, n. 3 del 2006), basata su una “Lezione
Lagrangiana” tenuta a Torino dallo stesso autore.
Per ulteriori approfondimenti sul caos deterministico rimandiamo a [F, capitolo
6] e anche all’articolo “Sistemi dinamici e caos deterministico” di Gian-Italo Bischi,
Modelli Matematici applicati all’Ecologia (3.12.06)
76
Capitolo 6
Rosa Carini, Laura Gardini e Paolo Tenti, pubblicato sul n. 47 di Lettera Matematica Pristem e disponibile alla pagine web
http://matematica.uni-bocconi.it/caos/home.htm
§ 6.3 Stabilità per punti di equilibrio di sistemi di equazioni
differenziali
Diamo ora una definizione più precisa di stabilità di un punto di equilibrio isolato.
Un punto di equilibrio (x0 , y0 ) è detto metastabile se le orbite che partono vicino ad esso, rimangono sempre vicino, ossia, se per ogni suo intorno
I, si può trovare un altro intorno I 0 , contenuto in
I e tale che ogni soluzione (x(t), y(t)) che passa da
un punto di I 0 , rimane in I, per ogni t successivo.
Alcuni autori chiamano stabili i punti che chiamiamo qui metastabili.
Un punto di equilibrio metastabile (x0 , y0 ) è
detto stabile (o secondo alcuni autori asintoticamente stabile) se per ogni suo intorno I, si
può trovare un altro intorno I 0 , contenuto in I e
tale che ogni soluzione (x(t), y(t)) che passa da un
punto di I 0 , tende al punto (x0 , y0 ), cioè
lim(x(t), y(t)) = (x0 , y0 ) .
t→0
Piano delle fasi e stabilità per sistemi lineari omogenei
Abbiamo visto come le soluzioni di un sistema di equazioni lineare omogeneo dipendano dagli autovalori della matrice A del sistema. In effetti i diversi tipi di piani
delle fasi e relative orbite dipendeno solo dagli autovalori.
S. Console – M. Roggero
Stabilità
ti
77
Osserviamo innanzitutto che per un sistema lineare omogeneo a coefficienti costan 0
x = ax + by
y 0 = cx + dy
l’origine è un punto stazionario.
Inoltre, se il determinante della matrice A è non nullo, cioè ad − bc 6= 0, l’origine
è l’unico punto stazionario.
Supporremo per semplicità che det A 6= 0, il che significa che gli autovalori di A
non sono nulli.
Per studiare qualitativamente il piano delle fasi possiamo sempre trasformare
le coordinate in modo che le direzioni degli autovettori (o almeno uno di essi)
coincidano con quelle degli assi. In tale riferimento, il sistema assume la forma
canonica
0
x = λ1 x
y 0 = λ2 y
se gli autovalori sono reali (distinti, o coincidenti, purchè esistano due autovettori
indipendenti), oppure un sistema lineare omogeneo a coefficienti costanti
0
x = αx − βy
y 0 = βx + αy
se gli autovalori sono complessi α ± iβ.
Esaminiamo il piano delle fasi e la stabilità dell’unico punto stazionario, distinguendo i vari casi.
Modelli Matematici applicati all’Ecologia (3.12.06)
78
Capitolo 6
ˆ A ha autovalori reali distinti (e non nulli). Distinguiamo 3 casi a seconda del
loro segno
1. Gli autovalori sono entrambi negativi :
abbiamo un nodo stabile e il piano delle fasi e relative orbite è come in
figura
2. Gli autovalori sono entrambi positivi :
abbiamo un nodo instabile. Come mostrato in figura, le orbite sono le
stesse dei caso precedente, ma “escono dall’origine, invece che entrarci”.
S. Console – M. Roggero
Stabilità
79
3. Gli autovalori hanno segno opposto:
abbiamo un punto di sella, che è instabile.
ˆ A ha autovalori coincidenti (e non nulli)
Una possibilità è che ci siano 2 autovettori indipendenti. In tal caso, nella base
formata dagli autovettori, la matrice si scrive come una matrice diagonale del
tipo
λ 0
A=
.
0 λ
In questo caso abbiamo un nodo stellato, che è stabile se l’autovalore λ è
negativo (come in figura), instabile se l’autovalore λ è positivo.
Se invece, l’autospazio relativo all’unico autovalore è di dimensione 1, in una
opportuna base (di cui un elemento è un autovettore, la matrice del sistema si
Modelli Matematici applicati all’Ecologia (3.12.06)
80
Capitolo 6
scrive nella forma
A=
λ 0
γ λ
.
Abbiamo allora un nodo degenere, che è stabile se l’autovalore λ è negativo
(come in figura), instabile se l’autovalore λ è positivo.
ˆ A ha autovalori complessi coniugati α ± iβ.
Se α = 0: abbiamo un centro o vortice che è metastabile ma non stabile.
In figura abbiamo un caso “generale” in cui le orbite sono ellissi (in figura
tuttavia gli assi delle ellissi sono gli assi coordinati: in generale posso avere
“ellissi ruotate”)
S. Console – M. Roggero
Stabilità
81
Un caso particolare si ha se le orbite sono circonferenze concentriche. Questo avviene ad
esempio se la matrice del sistema ha la forma
0 β
A=
.
−β 0
Se α 6= 0: abbiamo allora un fuoco che è stabile se α < 0 (caso della figura),
instabile se α > 0.
In conclusione abbiamo il seguente
Teorema 6.9. Sia dato un sistema lineare omogeneo a coefficienti costanti
0
x = ax + by
y 0 = cx + dy
Modelli Matematici applicati all’Ecologia (3.12.06)
82
Capitolo 6
e supponiamo che il determinante della matrice A sia non nullo, cioè ad − bc 6= 0.
Siano λ1 e λ2 gli autovalori di A, che sono non nulli. In generale, λ1 e λ2 sono
numeri complessi ed indichiamo con Re (λ1 ) e Re (λ2 ) le loro parti reali.
Allora l’origine è l’unico punto stazionario che è
ˆ stabile se Re (λ1 ) < 0 e Re (λ2 ) < 0;
ˆ instabile se Re (λ1 ) > 0 e/o Re (λ2 ) > 0;
ˆ metastabile se Re (λ1 ) = Re (λ2 ) = 0.
Stabilità per sistemi generali
Consideriamo il sistema generale
x0 = f (x, y)
y 0 = g(x, y)
per cui il punto P0 = (x0 , y0 ) sia un punto stazionario. Sotto opportune ipotesi di
regolarità per le funzioni f (x, y) e g(x, y) possiamo approssimare il sistema con un
sistema lineare. L’idea è di fare lo sviluppo di Taylor delle funzioni f (x, y) e g(x, y)
basato nel punto P0 = (x0 , y0 ) (cfr. Appendice). Per semplicità possiamo supporre
P0 = (x0 , y0 ) = (0, 0) e siccome P0 è un punto di equilibrio, f (x0 , y0 ) = g(x0 , y0 ) = 0.
Allora gli sviluppi di Taylor arrestati al prim’ordine di f e g sono dati da
f (x, y) = fx0 (0, 0)x + fy0 (0, 0)y + ω(x, y)
g(x, y) = gx0 (0, 0)x + gy0 (0, 0)y + ω 0 (x, y)
p
con ω(x, y) → 0 e ω 0 (x, y) → 0 per r = x2 + y 2 → 0.
Dunque possiamo approssimare il sistema generale con il sistema linearizzato
0
x = fx0 (0, 0)x + fy0 (0, 0)y
y 0 = gx0 (0, 0)x + gy0 (0, 0)y
che è un sistema lineare omogeneo con (genericamente) unico punto di equilibrio
l’origine O(0, 0).
La matrice del sistema linearizzato, è dunque la matrice Jacobiana
0
fx fy0
gx0 gy0
Ci domandiamo quando il sistema linearizzato ha le stesse proprietà del sistema
originario in merito alla stabilità del punto di equilibrio O(0, 0)?
Se ci ricordiamo (Esempio 5.5), abbiamo osservato che una perturbazione di un
centro (metastabile) produce una spirale che può essere sia stabile, sia instabile.
S. Console – M. Roggero
Stabilità
83
Dunque è lecito pensare che se il sistema linearizzato ha un centro nulla si possa
dire sul sistema originario non lineare. Negli altri casi si può invece usare il sistema
linearizzato per trarre conclusioni sulla stabilità del punto di equilibrio del sistema
non lineare.
Teorema 6.10. Supponiamo che il sistema non lineare
0
x = f (x, y)
y 0 = g(x, y)
abbia P0 = (x0 , y0 ) come punto stazionario isolato. Sotto opportune ipotesi di
regolarità per le funzioni f (x, y) e g(x, y), considerato il sistema linearizzato
0
x = fx0 (0, 0)x + fy0 (0, 0)y
y 0 = gx0 (0, 0)x + gy0 (0, 0)y
(traslando P0 nell’origine O), se il sistema linearizzato ha in O un fuoco, un nodo
o un punto a sella, tale punto è dello stesso tipo per il sistema non lineare e il suo
carattere di stabilità è il medesimo.
§ 6.4
Esercizi risolti
6.1 Determinare i punti di equilibrio del sistema
 0
x = −x + eαy − 1
y 0 = ex−y − 1
e studiarne la stabilità.
Soluzione: I punti di equilibrio sono le soluzioni del sistema

−x + eαy − 1 = 0
.
ex−y − 1 = 0
Dalla seconda equazione si trova
ex−y = 1 ,
x−y =0,
x=y
e sostituendo nella prima
−y − 1 + eαy = 0 ,
y + 1 = eαy ,
y=0.
Dunque l’unico punto di equilibrio è l’origine O(0, 0).
Il sistema linearizzato nell’origine è
„ 0 « „
«„
«
x
−1 α
x
=
.
y0
1
−1
y
Gli autovalori della matrice del sistema sono
√
1. λ = −1 ± α se α ≥ 0. In questo caso, se 0 < α < 1 gli autovalori sono non nulli e negativi
e dunque l’origine è un nodo stabile. Se α > 1 l’origine è una sella, mentre se α = 1, un
autovalore è nullo e non si può applicare il teorema di linearizzazione (Teorema 6.10).
√
2. λ = −1 ± i −α se α < 0 e l’origine è un fuoco stabile.
Modelli Matematici applicati all’Ecologia (3.12.06)
84
§ 6.5
Capitolo 6
Altri esercizi
6.2 Una popolazione di pesci in un lago segue un modello di crescita logistica, ma c’è un prelievo
costante dovuto alla pesca. Il tasso di crescita iniziale (“malthusiano”) è k = 0.2 per mese, la
capacità portante è di 40 (migliaia) di pesci e la pesca è di 1.5 (migliaia al mese). Scrivere il
modello di crescita, trovando e classificando gli equilibri. La popolazione di pesci si estingue? Qual
è la popolazione nel lungo termine?
6.3 Determinare e classificare i punti di equilibrio dei seguenti modelli continui:
x0 = x2 (3 − x) ,
x0 = 2x(1 − x) − 21 x ,
x0 = (4 − x)(2 − x)2 .
6.4 Determinare e classificare i punti di equilibrio dei seguenti modelli continui in dipendenza dal
parametro h; disegnare il diagramma di biforcazione:
x0 = hx − x2 ,
0
x = (1 − x)(x2 − h) .
6.5 Determinare i punti di equilibrio del sistema
 0
x = (2 − x − y)x
y 0 = (3 − 2x − y)y
e studiarne la stabilità.
S. Console – M. Roggero
Capitolo 7
Modelli di popolazioni conviventi
Una popolazione si può considerare isolata se non è possibile individuare, all’interno dell’ecosistema, un partner privilegiato, una popolazione cioè la cui presenza
o assenza determini direttamente una mutazione sulla composizione della prima
popolazione. Due o più popolazioni la cui crescita sia interdipendente si definiscono
conviventi. Sono popolazioni conviventi ad esempio una popolazione di predatori
con la popolazione della preda preferita, una popolazione di impollinatori ed una
popolazione di vegetali dal nettare particolare etc...
In questo capitolo studiamo quindi il caso di popolazioni conviventi. Si suppone
che le popolazioni considerate siano esenti da fenomeni migratori, che l’ambiente
in cui vivono non subisca mutazioni, che le popolazioni siano omogenee dal punto
di vista della natalità, riproduzione e mortalità (in altre parole che non vi siano
differenze tra individui giovani o vecchi) e che ogni popolazione risponda istantaneamente con una variazione a qualsiasi variazione dell’altra (anche se si verifica che in
realtà la variazione non è quasi mai istantanea ma spesso ritardata: una mancanza
di prede nel periodo estivo può portare ad una diminuzione dei predatori nel periodo
invernale, dove la malnutrizione fa si che gli individui, indeboliti, siano soggetti ad
una mortalità maggiore).
Le variazioni di grandezza di una popolazione sono dipendenti dal suo interagire con
tutte le specie presenti nel suo habitat.
Siccome sarebbe molto laborioso analizzare tutte queste interazioni insieme, per
esigenze di modellizzazione si attua una semplificazione e si fa riferimento all’interazione che avviene con una sola popolazione.
Chiamate X1 e X2 le popolazioni conviventi, denotiamo x1 = x1 (t) e x2 = x2 (t)
la loro numerosità al tempo t. Dalle ipotesi segue che le grandezze x1 e x2 devono
soddisfare un sistema di equazioni differenziali del tipo

dx1



 dt = x1 f1 (x1 , x2 )



 dx2 = x2 f2 (x1 , x2 )
dt
85
(3)
86
Capitolo 7
Si può notare che questa volta il tasso di accrescimento dipende dalla numerosità
di entrambe le popolazioni.
A rigor di logica, se si considera un biosistema in cui sia presente una popolazione
di carnivori, sarà da prendere in considerazione anche una popolazione di erbivori e,
di conseguenza, una popolazione di vegetali necessari al loro sostentamento. Per esigenze di semplicità, tuttavia, spesso si considera solo l’interazione carnivoro-erbivoro
tralasciando quella erbivoro-vegetale ipotizzando, cioè, una disponibilità infinita della popolazione vegetale, o supponendo che questa abbia un effetto diretto e limitativo
della crescita della popolazione erbivora.
Considerando gli effetti demografici determinati dalla reciproca influenza delle due
popolazioni si possono avere vari tipi di interazioni:
Mutualistiche (o cooperative) Si ha tra due popolazioni che risultano mutuamente avvantaggiate dall’interazione stessa. In questo caso si dovrà avere in
(3)
∂f1
>0
∂x2
e
∂f2
>0
∂x1
in quanto ad una crescita di X2 deve corrispondere una crescita di X1 e viceversa.
Amensalistiche Si ha tra due popolazioni delle quali una (sia X1 ) risulta penalizzata dall’interazione, mentre l’altra non ne risente. In questo caso si dovrà
avere in (3)
∂f1
<0
∂x2
e
∂f2
=0
∂x1
in quanto la presenza di X2 provoca una diminuzione di x1 mentre x2 non è
influenzata da X1 .
Commensalistiche Si ha tra due popolazioni delle quali una (sia X1 ) risulta avvantaggiata dall’interazione, mentre l’altra non ne risente. In questo caso si
dovrà avere in (3)
∂f1
>0
∂x2
e
∂f2
=0
∂x1
in quanto la presenza di X2 provoca un accrescimento di x1 mentre x2 non è
influenzato da X1 .
S. Console – M. Roggero
Modelli di popolazioni conviventi
87
Competitive Si ha tra due popolazioni che risultano mutuamente penalizzate dall’interazione stessa. In questo caso si dovrà avere in (3)
∂f1
<0
∂x2
e
∂f2
<0
∂x1
in quanto ad una crescita di X2 deve corrispondere una decrescita di X1 e
viceversa.
Di predazione Si ha tra due popolazioni delle quali una (sia X1 ), detta predatore,
trae vantaggio dall’interazione, mentre l’altra, detta preda, risulta penalizzata
dalla stessa. In questo caso si dovrà avere in (3)
∂f1
>0
∂x2
e
∂f2
<0
∂x1
in quanto ad una crescita di X2 deve corrispondere una crescita di X1 , mentre
ad una crescita di X1 deve corrispondere una decrescita di X2 .
Nel caso in cui l’interazione delle popolazioni non produca effetti sulla loro
crescita, si ha ovviamente
∂f1
∂f2
=
= 0.
∂x2
∂x1
§ 7.1
Modelli quadratici
Partendo dal sistema generale (3), Alfred Lotka e Vito Volterra (separatamente ma
giungendo alle medesime conclusioni) ne hanno studiato alcuni casi particolari: si
tratta di modelli quadratici associati al sistema di equazioni differenziali

dx1



 dt = x1 (a1 + a11 x1 + a12 x2 )
(4)



 dx2 = x2 (a2 + a21 x1 + a22 x2 )
dt
dove a1 , a2 , a11 , a12 , a21 , a22 sono costanti.
Questo sistema, di più semplice trattazione, si rivela strumento molto efficace nello
studio delle variazioni demografiche.
Modelli Matematici applicati all’Ecologia (30.07.06)
88
Capitolo 7
dx1
= 0 sulla coppia di rette x1 = 0 e
dt
dx2
a1 + a11 x1 + a12 x2 = 0 e analogamente
= 0 sulla coppia di rette x2 = 0 e
dt
a2 + a21 x1 + a22 x2 = 0. Tali rette sono dette isocline, perchè su ciascuna delle
coppie il campo ha la stessa direzione parallela agli assi.
Notiamo che dalla (4) si ricava che
I punti di intersezione delle isocline (prese a due a due, una per coppia di rette)
danno i punti di equilibrio
O = (0, 0) ,
E3 =
E1 =
a1
,0
a11
,
E2 =
a1 a22 − a2 a12 a2 a11 − a1 a21
,
a11 a22 − a12 a21 a11 a22 − a12 a21
S. Console – M. Roggero
a2
0,
a22
.
,
Modelli di popolazioni conviventi
§ 7.2
89
Predazione: il modello Lotka-Volterra
Il modello Lotka-Volterra propriamente detto è quello relativo al caso di interazione
preda-predatore.
Dall’analisi precedente, relativa al sistema (3), si deduce che, affinchè il sistema (4)
rappresenti effettivamente la situazione in esame, considerate X2 come popolazione
predatrice e X1 come preda, deve essere:
∂f1
∂(a1 + a11 x1 + a12 x2 )
=
<0
∂x2
∂x2
⇒ a12 < 0
∂f2
∂(a2 + a21 x1 + a22 x2 )
=
>0
∂x1
∂x1
⇒ a21 > 0
e
Le ipotesi che stanno alla base del modello in esame sono:
1. In assenza della popolazione X2 la popolazione X1 cresce esponenzialmente.
Se x2 = 0 deve quindi essere
dx1
= a1 x1
dt
con a1 > 0.
2. In assenza della popolazione X1 la popolazione X2 decresce esponenzialmente.
Se x1 = 0 deve quindi essere
dx2
= −a2 x2
dt
con a2 > 0.
3. In presenza della popolazione X2 il tasso di crescita di X1 deve diminuire
proporzionalmente ad x2 .
Possiamo quindi scrivere
dx1
= (a1 − a12 x2 )x1 .
dt
4. In presenza della popolazione X1 il tasso di crescita di X2 deve aumentare
proporzionalmente ad x1 .
Possiamo quindi scrivere
dx2
= (−a2 + a21 x1 )x2 .
dt
Modelli Matematici applicati all’Ecologia (30.07.06)
90
Capitolo 7
Le prime due ipotesi si basano sulla semplificazione che la dieta della popolazione
X2 sia costituita esclusivamente da individui della popolazione X1 , senza la quale
non potrebbe sopravvivere.
Per convincerci della bontà dell’equazione scritta nella terza, possiamo osservare che
il controllo effettuato dal predatore sul tasso di crescita della preda è proporzionale
al numero di incontri che avvengono tra prede e predatori ed al numero di questi
che si concludono con l’eliminazione della preda stessa. E’ ragionevole supporre che
il numero di incontri sia proporzionale al numero di individui delle due popolazioni,
sia αx1 x2 , e che il numero di uccisioni sia proporzionale al numero di incontri, sia
βαx1 x2 . Ponendo βα = a12 si ottiene la terza ipotesi.
Con considerazioni analoghe si ottiene l’equazione scritta nell’ultima ipotesi.
Riassumendo, l’interazione tra preda e predatore è modellata dal sistema di equazioni
differenziali

dx1



 dt = (a1 − a12 x2 )x1
(5)



 dx2 = (−a2 + a21 x1 )x2
dt
dette equazioni di Lotka-Volterra preda-predatore.
Studio delle equazioni di Lotka-Volterra
Prima di cominciare lo studio delle equazioni da un punto di vista strettamente
matematico, può essere utile cercare di capire che tipo di risultato ci si possa aspettare.
In questo modello possiamo osservare come i predatori possano proliferare quando
vi sono molte prede, ma, a lungo andare, sovrasfruttano la loro fonte di cibo e vanno
incontro ad un declino. Siccome diminuisce la popolazione dei predatori, le prede
tenderanno ad aumentare di numero, gettando le basi per un nuovo aumento della
popolazione dei predatori in un ciclo di crescite e declini 1 con andamento simile a
quello di figura.
1
Queste sono da ritenersi valide, ovviamente, se il numero di individui delle due popolazioni non
scende al di sotto del numero minimo vitale: in quel caso la popolaizone in oggetto andrà incontro all’estinzione indipendentemente dallo sviluppo dell’altra, che subirà tuttavia le conseguenze
dell’estinzione del proprio partner privilegiato.
S. Console – M. Roggero
Modelli di popolazioni conviventi
91
Dalla seguente figura si vede come questo andamento oscillante delle densità di
popolazione si rifletta nelle corrisponti orbite nel piano delle fasi.
Iniziamo lo studio vero e proprio con l’analisi dei punti di equilibrio.
Si ha equilibrio quando il numero di individui di ciascuna popolazione si mantiene
Modelli Matematici applicati all’Ecologia (30.07.06)
92
Capitolo 7
costante, laddove cioè il sistema (5) diventa

dx1



 dt = 0
(6)



 dx2 = 0
dt
cioè
(a1 − a12 x2 )x1 = 0
(−a2 + a21 x1 )x2 = 0
(7)
che ha come soluzioni
x1 = x2 = 0
e
a2
a1
x1 =
, x2 =
a21
a12
Ci sono quindi due posizioni di equilibrio delle quali la prima rappresenta l’estinzione
di entrambe le specie e la seconda rappresenta un punto fisso nel quale entrambe le
popolazioni mantengono il loro valore non nullo.
Il valore assunto da x1 e x2 dipende ovviamente dalle condizioni iniziali.
Studiamo ora la stabilità di questi punti di equilibrio. La matrice Jacobiana del
modello preda predatore è:
a1 − a12 x2
−a12 x1
J=
a21 x2
a21 x1 − a2
Valutata in (0, 0) diviene
J(0, 0) =
a1
0
0 −a2
Gli autovalori di questa matrice sono
λ1 = a1
λ2 = −a2
S. Console – M. Roggero
Modelli di popolazioni conviventi
93
Abbiamo già sottolineato che, affinchè il modello sia sensato, a1 e a2 devono
essere entrambi maggiori di zero. Da ciò discende che i due autovalori abbiano sempre segni differenti e quindi, come ci aspettavamo, l’origine è un punto di equilibrio
instabile, più precisamente un punto di sella.
Se l’origine fosse stata un punto di equilibrio stabile, infatti, le popolazioni iniziali
sarebbero state attratte da questo punto e quindi l’evoluzione del sistema avrebbe
portato all’estinzione di entrambe le specie per molte coppie di valori iniziali.
Passiamo ora allo studio del secondo punto di equilibrio.
Sostituendo le sue coordinate nella matrice Jacobiana si ottiene:
J
a2 a1
,
a21 a12
=
0
a1 a21
a12
− aa1221a2
0
Gli autovalori di questa matrice sono
√
λ1 = i a1 a2 ,
√
λ2 = −i a1 a2
Poichè gli autovalori sono entrambi complessi con parte reale nulla, questo punto di
equilibrio è un centro. Questo è proprio uno dei casi eccezionali in cui non possiamo
trarre alcuna conclusione dalla stabilità del sistema linearizzato (cf. Teorema 6.10):
per il sistema del modello Lotka-Volterra il punto può a priori essere una spirale
stabili o instabile, oltre che un centro. Si può però vedere con un calcolo diretto
(cui rimandiamo al testo di Capelo già citato) che si ha effettivamente un centro.
Riportiamo un’immagine ottenuta al computer di orbite vicine al punto di equilibrio:
Dunque se siamo vicini al punto di equilibrio, la
grandezza delle due popolazioni varia ciclicamente
oscillando intorno a questo
punto.
Il ritratto di fase di questo sistema di equazioni differenziali sarà quindi del tipo
Modelli Matematici applicati all’Ecologia (30.07.06)
94
Capitolo 7
È interessante dare un’interpretazione delle orbite in figura relativamente alla
grandezza delle popolazioni.
Le orbite che coincidono con i semiassi rappresentano le prime due ipotesi in quanto
ad un’assenza del predatore (x2 = 0) segue una crescita delle prede tendente a
infinito e ad un’assenza di prede (x1 = 0) corrisponde il collasso della popolazione
di predatori.
Osserviamo ora le orbite iperboliche. Se ad un dato istante il rapporto tra numero
di predatori e numero di prede è troppo sbilanciato in favore dei predatori, la caccia
risulta molto difficile e x2 scende al di sotto del numero minimo necessario alla
sopravvivenza della popolazione. A causa di questa diminuzione il numero di prede
aumenta, ma la popolazione predatrice è comunque destinata all’estinzione.
Infine rimangono le orbite chiuse. Consideriamo l’orbita O2 . Supponiamo che il
1
numero di predatori sia uguale alla condizione di equilibrio (cioè x2 = aa12
) e il
numero di prede sia invece maggiore, e chiamiamo questo punto P2 . Poichè vi è
abbondanza di prede, x2 tende ad aumentare, provocando una diminuzione di x1 .
2
Arrivati a M1 , punto di equilibrio per il numero di prede (cioè x1 = aa21
), anche
a1
x2 comincia a diminuire, per scarsità di cibo. Quando x2 = a12 (punto M2 ), x1
ha nuovamente la possibilità di aumentare (in quanto la probabilità di incontro
preda-predatore è diminuita) fino a giungere ad M3 (equilibrio per le prede), dove
l’abbondanza relativa di prede provoca un nuovo aumento di x2 fino a giungere in
P2 all’equilibrio.
S. Console – M. Roggero
Modelli di popolazioni conviventi
95
Un’applicazione del modello: effetto pesca
Un esempio classico di applicazione del modello Lotka-Volterra è il caso studiato
da Volterra stesso basandosi su dati sperimentali fornitigli dallo zoologo Umberto
D’Ancona.
In quel periodo la pesca era effettuata con reti e il D’ancona registrava la percentuale
di pesci della classe dei selaci pescati, in quanto questi ultimi sono predatori delle
specie cui era rivolta in realtà la caccia. Per utilizzare il modello si deve modellare
l’effetto della pesca sul biosistema costituito dai selaci e dagli altri pesci, in questo
contesto predatori e prede.
Consideriamo due tipi di prelevamento:
prelevamento una tantum in un dato istante viene prelevata una quantità q1 di
individui della popolazione X1 e q2 della popolazione X2
prelevamento continuo in ogni istante viene prelevata una frazione p1 x1 della
popolazione X1 ed una frazione p2 x2 della popolazione X2
Osserviamo che nel secondo caso può essere opportuno scrivere i coefficienti di
prelevamento nella forma p1 = λπ1 e p2 = λπ2 , dove λ rappresenta l’intensità del
prelevamento e πi le modalità di prelevamento.
Nel caso del prelevamento continuo, che è quello realmente utilizzato nell’analisi di
D’Ancona e Volterra, il sistema (5) diventa

dx1



 dt = (a1 − p1 − a12 x2 )x1
(8)



 dx2 = (−a2 − p2 + a21 x1 )x2
dt
La sua Jacobiana è
J=
a1 − p1 − a12 x2
−a12 x1
a21 x2
−a2 − p2 + a21 x1
Per quanto riguarda la seconda equazione del sistema, la modifica non è sostanziale
in quanto −p2 introduce semplicemente un’ulteriore diminuzione del tasso di crescita
della popolazione predatrice, già di per se negativo. Invece la modifica nella prima
equazione può portare ad una variazione anche qualitativa della dinamica in oggetto. Al variare dell’intensità di prelevamento, quindi, possono presentarsi tre diverse
situazioni.
Modelli Matematici applicati all’Ecologia (30.07.06)
96
Capitolo 7
Caso I: a1 − p1 < 0
Per trovare gli equilibri bisogna risolvere il sistema
(a1 − p1 − a12 x2 )x1 = 0
(−a2 − p2 + a21 x1 )x2 = 0
che, posto k1 = |a1 − p1 | e k2 = | − a2 − p2 |, diventa
(−k1 − a12 x2 )x1 = 0
(−k2 + a21 x1 )x2 = 0
Le soluzioni sono
x1 = x2 = 0
e
k2
k1
x1 =
, x2 = −
a21
a12
La seconda coppia di valori è
da scartare in quanto non è possibile avere una popolazione formata da un numero negativo di
individui.
La Jacobiana nell’unico punto di
equilibrio è
a1 − p1
0
J(0, 0) =
0
−a2 − p2
che ha autovalori
λ1 = a1 − p1 < 0 ,
λ2 = −(a2 + p2 ) < 0 .
Dunque l’origine è un nodo stabile, pertanto sia la popolazione preda che la popolazione predatore tendono all’estinzione qualunque siano le condizioni iniziali.
Caso II: a1 − p1 = 0
In questo caso il sistema (8) diventa

dx1



 dt = −a12 x1 x2



 dx2 = (−a2 − p2 + a21 x1 )x2
dt
S. Console – M. Roggero
Modelli di popolazioni conviventi
97
ed ha come punti di equilibrio tutti i punti del semiasse x1 ≥ 0, che non sono
quindi isolati.
Le orbite del sistema saranno allora le soluzioni dell’equazione differenziale
dx2
−a2 − p2 + a21 x1
=
dx1
−a12 x1
avendo posto la condizione x1 > 0. Con la condizione iniziale x20 = x2 x10
x2 (x1 ) =
a2 + p2
x1
a21
log (
)−
(x1 − x10 ) + x20
a12
x10
a12
è un’orbita del nostro sistema. Possiamo osservare che tutte le soluzioni hanno
2
un massimo per x1 = a2a+p
.
21
Diamo ora un’interpretazione di questo grafico.
Se il sistema si trova nella condizione x2 = 0, cioè non vi sono predatori, la popolazione X1 si mantiene costante, come d’altra parte era lecito aspettarsi dato che
a1 − p1 rappresenta il tasso di crescita della popolazione isolata.
Se invece il sistema si trova in un punto P non appartenente all’asse delle ascisse,
data la presenza di prede la popolazione X2 tenderà ad ingrandirsi fino a raggiungere un massimo per poi diminuire fino all’estinzione. Nello studio dell’equazione
differenziale avevamo posto x1 > 0, ma è ovvio che nella condizione x1 = 0 (e solo
sotto questa ipotesi) si giungerebbe all’estinzione di entrambe le specie.
Modelli Matematici applicati all’Ecologia (30.07.06)
98
Capitolo 7
Caso III: a1 − p1 > 0
Ponendo k1 = a1 − p1 e k2 = −a2 − p2 , si ottiene un sistema analogo al sistema
(5) già analizzato.
Autolimitazione della preda
Mutua interferenza
Risposta funzionale del predatore
Nel modello preda-predatore di Lotka-Volterra “classico” (LVPP)
oppure con autolimitazione della preda
il termine di predazione è supposto proporzionale al prodotto delle densità x1 della
preda e x2 del predatore.
Questo non è realistico in molte situazioni, in quanto presuppone che la capacità
di predazione sia limitata solo dalla disponibilità di prede e quindi che il predatore
possa cacciare e consumare prede senza limiti.
Viene detta risposta funzionale del predatore quantità di preda utilizzata
H in funzione delle prede totali P disponibili (in un dato tempo).
Ci sono 3 tipi di risposte funzionali, illustrate dai seguenti grafici:
S. Console – M. Roggero
Modelli di popolazioni conviventi
99
ˆ Risposta funzionale di tipo 1:
La funzione è lineare all’inizio (quindi il termine moltiplicativo in LVPP) poi ho una saturazione. Questo tipo di risposta non è molto
comune ed è limitata agli animali filtratori (ad
esempio la Daphnia magna rispetto al lievito
Saccharomyces cerevisiae)
ˆ Risposta funzionale di tipo 2: La velocità di estrazione aumenta al crescere
delle risorse ma con andamento gradualmente decelerato. Siano
T : tempo totale dedicato alla nutrizione
Tr : tempo dedicato alla ricerca della preda
Tm : tempo dedicato alla manipolazione della preda
tm : tempo di manipolazione per unità di preda
q : volume esplorato dal predatore in unità di tempo
pa : probabilità che il predatore catturi la preda
Modelli Matematici applicati all’Ecologia (30.07.06)
100
Capitolo 7
Il numero di prede incontrate in Tr è qx1 Tr , quindi la quantità di prede catturate in Tr è C = pa qx1 Tr . Per manipolarle occorre un tempo Tm = tm C =
tm pa qx1 Tr . D’altra parte T = Tr +Tm = Tr +tm pa qx1 Tr e dunque il numero di
prede assimilate nel tempo T è H = C/T = pa qx1 T/ (Tr + tm pa qx1 Tr ). Allora
H=
pa qx1
.
1 + tm pa qx1
Dato che la preda trattata nel tempo T è P = pa qx1 ,
H=
P
.
1 + tm P
Posto p = pa q e t = tm pa qx1 , il modello LVPP
con autolimitazione della preda diventa

px1 x2
dx1


= a1 (1 − x1 /K) −

 dt
1 + tx1



 dx2 = −a2 x2 + e px1 x2
dt
1 + tx1
L’analisi di queste equazioni non è semplice e
ho un diverso comporatamento al variare dei
parametri (v.ad es. [ML, pag. 57]).
ˆ Risposta funzionale di tipo 3: Nella risposta di tipo 2 si è supposto che i
parametri di predazione e in particolare pa e q fossero costanti.
In realtà all’aumentare della densità “l’interesse del predatore verso la preda” aumenta,
creando un effetto valanga per densità crescenti della preda, ma poi subentra un effetto
di saturazione come nella risposta di tipo 2.
Allora ho una curva sigmoide
H=
ap
1 + bP + tm P 2
S. Console – M. Roggero
Modelli di popolazioni conviventi
§ 7.3
101
Competizione interspecifica
Per popolazioni in competizione, la coesistenza può essere:
ˆ stabile oppure
ˆ di esclusione competitiva che a sua volta può essere
– condizionale cioè dipendente dalle condizioni iniziali;
– assoluta.
Modello Lotka-Volterra quadratico per la competizione di due specie
Ipotesi:
ˆ la velocità di accrescimento di una specie è funzione decrescente della sua
densità x1 della densità dell’altra specie x2 ;
ˆ idem per l’altra specie x2 ;
ˆ la crescita sia logistica.
Riscrivendo l’equazione logistica per singola specie come
dx
K −x
=x·r
,
dt
K
ed esprimendo la competizione intraspecifica di ciascuna specie con a11 e a22 e la
competizione interspecifica tra ciascuna specie con a12 e a21 si ha

K1 − a11 x1 − a12 x2
dx1


= x1 · r1
,


 dt
K1


dx2
K2 − a21 x1 − a22 x2


= x2 · r2
.

dt
K2
Se normalizziamo i coefficienti della competizione interspecifica alla competizione
intraspecifica

K1 − x1 − a12 x2
dx1


= x1 · r1
,


 dt
K1


dx2
K2 − a21 x1 − x2


=
x
·
r
.

2
2
dt
K2
Le isocline sono:
Modelli Matematici applicati all’Ecologia (30.07.06)
102
Capitolo 7
Sono possibili 3 casi:
1. L’isoclina zero di una delle popolazioni giace completamente nel semipiano dell’altra.
In questo caso sono a sua volta possibili vari casi:
1.a: la popolazione 2 è svantaggiata sulla
popolazione 1
K1
> K2
a12
K2
< K1
a21
Il sistema evolve verso la popolazione avvantaggiata 1: in questo caso si parla di
esclusione competitiva assoluta
S. Console – M. Roggero
Modelli di popolazioni conviventi
103
1.b: la popolazione 1 è svantaggiata
sulla popolazione 2
K1
< K2
a12
K2
> K1
a21
Il sistema evolve verso la popolazione
avvantaggiata 2: ovviamente questo
è il caso simmetrico di esclusione
competitiva assoluta
2. Le isocline si intersecano in un punto E che giace sotto la congiungente delle
capacità portanti K1 e K2 .
Per entrambe le popolazioni,
i vincoli per competizione interspecifica sono maggiori di
quelli per competizione intraspecifica.
Il sistema evolve verso K1
o K2 a seconda delle condizioni iniziali (la congiungente E con l’origine fa da
spartiacque).
In questo caso si ha una situazione di esclusione competitiva condizionale.
3. Le isocline si intersecano in un punto E che giace sopra la congiungente delle
capacità portanti K1 e K2 .
Modelli Matematici applicati all’Ecologia (30.07.06)
104
Capitolo 7
K1
> K2
a12
K2
> K1
a21
Per entrambe le popolazioni, i vincoli per competizione intraspecifica sono maggiori di quelli per
competizione interspecifica.
IIl sistema evolve verso E indipendentemente dalla condizione iniziale.
In questo caso si ha coesistenza
stabile.
§ 7.4
Cooperazione o mutualismo
Con “mutualismo” si intende il rapporto tra due individui di specie diverse che
traggono entrambi vantaggi dal loro vivere insieme.
ˆ entrambe le popolazioni sono dipendenti l’una dall’altra;
ˆ spesso sono associati degli organismi molto diversi tra loro (organismi che
hanno richieste molto differenti);
ˆ gli individui di ciascuna specie crescono e sopravvivono meglio e/o si riproducono con maggior successo quando sono in associazione;
ˆ spesso il vantaggio riguarda le risorse alimentari per almeno una specie, la
protezione o l’ambiente favorevole alla riproduzione per l’altra;
ˆ il mutualismo si è evoluto perchè i benefici che derivano a ciascun partner
superano i costi che l’associazione implica;
Tra gli effetti del mutualismo si ha
ˆ La difesa dai predatori: ad esempio questo succede per pesci pagliaccio (Amphiprion) ed attinie delle barriere coralline.
ˆ Aumenta la disponibilità di risorse: coltivazione di funghi basidiomiceti da
parte di formiche tagliatrici di foglie (Atta).
S. Console – M. Roggero
Modelli di popolazioni conviventi
105
Sorprendente associazione tra funghi Basidiomiceti e formiche della specie Atta cephalotes.
Le formiche raccolgono e sminuzzano foglie e petali di fiori sui quali crescono funghi particolari che si trovano esclusivamente all’interno dei formicai dove trovano un ambiente ideale
alla loro diffusione.
ˆ Consente agli animali di utilizzare lignina, cellulosa ed emicellulose: batteri
del rumine; batteri e protozoi dell’intestino delle termiti
ˆ Crea un ambiente protetto per uno dei mutualisti: questo capita in molti casi
precedenti, ed inoltre per Hydra viridis e Chlorella.
Modelli Matematici applicati all’Ecologia (30.07.06)
106
Capitolo 7
L’Hydra viridis possiede una colorazione verdognola e vive quasi in simbiosi alle alghe verdi del
genere Chlorella che si localizzano nei suoi tessuti.
ˆ Aumenta la competitività per i nutrienti: piante con noduli radicali (leguminoseRhizobium; attinomiceti-Alnus e Dryas).
Noduli radicali di Attinomicete in una pianta
di Ontano (Alnus): questi
noduli fissano l’azoto atmosferico.
ˆ Aumenta la capacità di assorbimento radicale (micorrizze).
S. Console – M. Roggero
Modelli di popolazioni conviventi
107
I funghi non sono soltanto organismi decompositori.
Almeno un quarto di essi (micorrizze)
instaurano relazioni mutualistiche con
le piante.
Sebbene il corpo fruttifero (carpoforo)
sia una struttura di dimensioni relativamente modeste, la restante parte sotterranea dl fungo (micelio) risulta molto
sviluppata.
Questo avvolge le radici delle piante,
determinandone un aumento della capacità di assorbimento. La pianta dal canto suo fornisce al fungo zuccheri e amido che
utilizza nel suo metabolismo.
ˆ Consente la riproduzione e la dispersione delle piante: impollinatori ed inseminatori.
ˆ Dà origine a organismi nuovi: i licheni. I licheni sono complesse associazioni
mutualistiche tra funghi e alghe. I funghi provvedono a fornire un sistema
meccanico di supporto (alle rocce, al legno, ecc.), e di protezione contro la
disidratazione per le alghe, oltre ad assorbire i nutrienti dal substrato. Le
alghe a loro volta utilizzano parte dei nutrienti forniti dal fungo e la CO2
da questo prodotta, nel processo di fotosintesi producono zuccheri e ossigeno
consumati a loro volta dal fungo.
Un’epifitica Myrmecodia ant-plant.
La base gonfia di questa pianta è cava, con
un’intricata galleria di camere in cui vivono
le formiche: le formiche difendono la pianta
e la pianta fornisce loro nutrienti.
Nello studio dei modelli di due popolazioni interagenti fra loro sono rari i testi
che parlano del mutualismo.
Questo è dovuto a due ragioni:
ˆ Ragione storica: il modello di Lotka-Volterra studia prevalentemente ed è stato
introdotto per i casi della competizione e di preda-predatore.
Modelli Matematici applicati all’Ecologia (30.07.06)
108
Capitolo 7
ˆ Ragione ecologica: si è osservato che il mutualismo è meno diffuso nelle regioni
temperate e boreali. Nelle regioni tropicali invece, con temperature più elevate
e foreste sempreverdi, le interazioni mutualistiche sono maggiori.
– Assenza di ant-plant mutualismi a N e S di 24◦ .
– Assenza di pipistrelli frugivori e nettarivori a N di 33◦ .
– Assenza di Orchid bee a N di 24◦ in America.
Cominciamo ora ad esaminare modelli matematici per lo studio del mutualismo.
Per due popolazioni la cooperazione può essere cosı̀ distinta (in prima approssimazione - vedremo in seguito una distinzione più fine):
ˆ Obbligatoria: entrambe le popolazioni si estinguono in assenza dell’altra.
ˆ Facoltativa: le popolazioni sono in grado di sopravvivere da sole, ma la
presenza di ciascuna favorisce la crescita dell’altra.
Modello quadratico per la cooperazione obbligatoria
Ipotesi:
ˆ La velocità di accrescimento di una specie è funzione decrescente della sua
densità: in assenza di x2 , x1 decresce in modo logistico:
dx1
= (−a1 − a11 x1 )x1 .
dt
ˆ Idem per la specie x2
dx2
= (−a2 − a22 x2 )x2 .
dt
ˆ La presenza di x2 fa crescere x1 proporzionalmente a x2 :
dx1
= (−a1 − a11 x1 + a12 x2 )x1 .
dt
ˆ Idem per l’altra specie
dx2
= (−a2 + a21 x1 − a22 x2 )x2 .
dt
In definitiva il modello corrisponde al sistema di equazioni:


 dx1 = (−a1 − a11 x1 + a12 x2 )x1 ,

 dt



 dx2 = (−a2 + a21 x1 − a22 x2 )x2 .
dt
S. Console – M. Roggero
Modelli di popolazioni conviventi
I punti di equilibrio sono
due:
– l’origine E1 , che è un nodo
stabile (che chiaramente corrisponde all’estinzione di entrambe le specie);
– il punto di incontro E2 delle
isoscline
−a1 − a11 x1 + a12 x2 = 0 e
−a2 + a21 x1 − a22 x2 = 0, che
è un punto di sella.
Modello quadratico per la cooperazione facoltativa
Ipotesi:
ˆ In assenza di x2 , x1 decresce in modo logistico:
dx1
= (a1 − a11 x1 )x1 .
dt
ˆ Idem per la specie x2
dx2
= (a2 − a22 x2 )x2 .
dt
ˆ La presenza di x2 fa crescere x1 proporzionalmente a x2 :
dx1
= (a1 − a11 x1 + a12 x2 )x1 .
dt
ˆ Idem per l’altra specie
dx2
= (a2 + a21 x1 − a22 x2 )x2 .
dt
In definitiva il modello corrisponde al sistema di equazioni:

dx1



 dt = (a1 − a11 x1 + a12 x2 )x1 ,



 dx2 = (a2 + a21 x1 − a22 x2 )x2 .
dt
Modelli Matematici applicati all’Ecologia (30.07.06)
109
110
Capitolo 7
Ci sono ora 3 casi per isocline e piano delle fasi, a seconda del segno di
D := det(aij ) = a11 a22 − a12 a21 .
Se D < 0,
le isocline oblique non si intersecano nel primo
quadrante:
l’origine è un nodo instabile;
i punti di incontro delle isocline “oblique” con
quella orizzontale, E1 e verticale E2 sono punti
di sella.
Se D > 0,
le isocline oblique si intersecano nel primo quadrante:
il loro punto di intersezione E3 è un nodo stabile;
l’origine è un nodo instabile;
i punti di incontro delle isocline “oblique” con
quella orizzontale, E1 e verticale E2 sono punti
di sella.
Se D = 0,
le isocline oblique sono parallele:
l’origine è un nodo instabile;
i punti di incontro delle isocline “oblique” con
quella orizzontale, E1 e verticale E2 sono punti
di sella.
S. Console – M. Roggero
Modelli di popolazioni conviventi
111
Dalle equazioni logistiche generiche di Lotka e Volterra:

dx1



 dt = (±a1 − a11 x1 + a12 x2 )x1 ,



 dx2 = (±a2 + a21 x1 − a22 x2 )x2 .
dt
in relazione ai segni di a1 e a2 si osserva che abbiamo in realtà 3 tipi di cooperazione
segni
++
−−
+−
tipo di mutualismo
facoltativo/facoltativo
obbligatorio/obbligatorio
facoltativo/obbligatorio
Limite dell’applicazione del modello quadratico per il mutualismo
Il modello quadratico è un modello ottimo per comprendere le dinamiche nelle interazioni preda-predatore e i meccanismi di competizione. Tuttavia è inadeguato per
capire le interazioni mutualistiche
Porta a soluzioni di crescita “esplosiva” di entrambe le specie con un insieme di
comuni benefici.
Solo nel caso “migliore” (cooperazione facoltativa con D > 0) non si prevede
crescita illimitata ne’ l’estinzione delle due specie.
Il modello non è dunque realistico e non è in accordo con i dati sperimentali.
Come deve essere un modello matematico che descriva la cooperazione tra due
specie?
Deve prevedere che ci sia la saturazione della crescita per almeno una delle
popolazioni
Come si può realizzare questo?
Sostituendo le funzioni del modello quadratico del tipo Lotka Volterra, della
dx
dx
forma
= (LINEARE) x con nuove funzioni
= (NON LINEARE) x proposte
dt
dt
nel modello di Kolmogorov
Possiamo descrivere il comportamento qualitativo dei sistemi mutualistici storcendo
almeno una delle isocline le quali corrispondono a funzioni del tipo:
f (x, y) = [a1 y/(1 + b1 x)] − c1 ,
g(x, y) = d2 (K − y) + [a2 x/(1 + b2 x)] .
Nel caso del mutulaismo facoltattivo-facoltativo le isoscline e il piano delle fasi nel
modello quadratico, si presenta nel seguente modo (in quest figure e le seguenti le
popolazioni sono indicate con N1 ed N2 invece che x1 e x2 ).
Modelli Matematici applicati all’Ecologia (30.07.06)
112
Capitolo 7
“Stocendo le isocline” si ha
S. Console – M. Roggero
Modelli di popolazioni conviventi
113
Nel caso della cooperazione obbligatoria-obbligatoria il cambiamento è il seguente
o rispettivamente
Modelli Matematici applicati all’Ecologia (30.07.06)
114
Capitolo 7
Nel caso della cooperazione obbligatoria-facoltativa, poniamo N1 essere la specie
obbligatoria cosı̀ che la isoclina incontri l’asse y a (0, a) ed N2 essere la specie
facoltativa cosı̀ che la isoclina incontri l’asse y in (0, M ).
Se a > M
oppure
S. Console – M. Roggero
Modelli di popolazioni conviventi
Se a < M
Modelli Matematici applicati all’Ecologia (30.07.06)
115
Appendice 1
Le successioni
Definizione 1.1. Si dice successione una qualsiasi funzione a : N → R. Spesso
per indicare una successione si usa la sequenza delle immagini:
a0 = a(0), a1 = a(1), . . . , an = a(n), . . . .
Talvolta si considera l’insieme N privato dello zero e quindi si considera come
primo termine di una successione quello con indice 1 ossia:
a1 = a(1), a2 = a(2), . . . , an = a(n), . . . .
Per questo motivo è opportuno fare attenzione al significato di espressioni come
“primo termine della successione ”oppure “i primi 5 termini della successione ”che
potrebbero risultare ambigui.
A priori qualsiasi sequenza di numeri costituisce una successione; potremmo ad
esempio costruirne una lanciando un dato e considerando come an il numero uscito
all’n-esimo lancio. In generale però saranno più semplici da trattare e spesso anche
più interessanti successioni i cui termini sono ottenibili mediante una qualche formula
matematica.
Esempio 1.2. La successione an = n2 è la successione dei quadrati dei numeri
naturali:
a0 = 0, a1 = 1, a2 = 4, a3 = 9, a4 = 16, . . . .
√
La successione bn = 3 n è la successione delle radici cubiche dei numeri naturali:
√
√
√
3
3
3
a0 = 0 , a1 = 1 , a2 = 2 , a3 = 3 , a4 = 4, . . . .
Come si vede dall’esempio precedente, due successioni diverse possono avere
alcuni termini uguali.
Un modo spesso usato per assegnare una successione è quello ricorsivo che
consiste nell’assegnare alcuni termini iniziali (il primo, oppure i primi due, oppure
. . . ) e una formula che permette di ottenere ogni termine successivo mediante quello
che lo precede (oppure quelli che lo precedono).
116
Le successioni
117
Esempio 1.3. La successione data da a0 = 3 e an = 2an−1 − 1 è:
a0 = 3 , a1 = 2 · 3 − 1 = 5 , a2 = 5 · 2 − 1 = 9 , a3 = 9 · 2 − 1 = 17, . . . .
Esempio 1.4. La successione di Fibonacci è data da a0 = 1, a1 = 1 e an =
an−2 + an−1 per ogni n ≥ 2:
a0 = 1 , a1 = 1 , a2 = 1+1 = 2 , a3 = 1+2 = 3 , a4 = 2+3 = 5 , a5 = 3+5 = 8 . . . .
A questo proposito si veda anche l’Appendice....
§ A. 1.1
Successioni aritmetiche
Definizione 1.5. Si dice successione (o progressione) aritmetica di termine
iniziale a0 e ragione d (con a0 , d ∈ R) la funzione a : N → R cosı̀ definita:
a(n) = an = a0 + nd.
Esplicitandone le immagini in sequenza, si scrive:
a0 , a1 = a0 + d , a2 = a0 + 2d , a3 = a0 + 3d , . . . .
La successione aritmetica an = a0 +nd si esprime in forma ricorsiva assegnando
il termine iniziale a0 e per ogni (n ≥ 1) la formula:
an = an−1 + d.
Abbiamo quindi un modo semplice di caratterizzare le successioni aritmetiche:
Una successione an è aritmetica ⇐⇒ la differenza tra due termini
consecutivi è costante ossia an − an−1 = d.
Esempio 1.6. La successione dei numeri pari 0, 2, 4, 6, . . . è la successione aritmetica di termine iniziale 0 e ragione 2, definita dalla legge an = 2n.
In forma ricorsiva possiamo assegnare questa successione come an = an−1 + 2,
(n ≥ 1), specificando che il termine iniziale a0 vale 0.
Esempio 1.7. La successione an = n2 dei quadrati dei numeri naturali non è una
successione aritmetica poiché ad esempio a4 − a3 = 16 − 9 = 7 mentre a2 − a1 =
4 − 1 = 3 6= 7.
Se sappiamo che una certa successione an è aritmetica, allora la conoscenza di
due qualsiasi dei suoi termini permette di determinare l’intera successione. Per ogni
coppia di numeri naturali distinti h, k si ha infatti:
ah − ak = (a0 + hd) − (a0 + kd) = (h − k)d
e quindi
d=
ah − ak
h−k
e
a0 = ah − hd.
Modelli Matematici applicati all’Ecologia (03.04.06)
118
Appendice 1
Esempio 1.8. La successione aritmetica an tale che a3 = 5 e a7 = 21 ha ragione
d = 21−5
7−3 = 4 e quindi a0 = a3 − 3d = −7 ossia è data da an = −7 + 4n.
Può essere utile ricordare la formula che dà la somma dei primi n termini di una
tale successione aritmetica:
a0 + a1 + · · · + an = (n + 1)a0 +
n2 + n
d
2
oppure, iniziando da n = 1:
a1 + · · · + an = na0 +
§ A. 1.2
n2 + n
d.
2
Successioni geometriche
Definizione 1.9. Si dice successione (o progressione) geometrica di termine
iniziale a0 e ragione q ( q ∈ R ) la funzione a : N → R cosı̀ definita:
a(n) = an = a0 q n .
Esplicitandone le immagini in sequenza, si scrive:
a0 , a1 = a0 q, a2 = a0 q 2 , a3 = a0 q 3 , . . . .
La successione geometrica an = a0 q n si esprime in forma ricorsiva ponendo an =
an−1 q, (n ≥ 1) e assegnando a0 come termine iniziale.
Abbiamo quindi un modo semplice di caratterizzare le successioni geometriche:
Una successione an è geometrica ⇐⇒ il rapporto tra due termini
an
consecutivi è costante ossia an−1
= q.
Esempio 1.10. La successione delle potenze di 2, ossia : 1, 2, 4, 8, 16, . . . è la
successione geometrica di termine iniziale 1 e ragione 2, definita dalla legge a(n) =
2n .
In forma ricorsiva possiamo scrivere an = 2an−1 , (n ≥ 1), specificando che il
termine iniziale a0 vale 1.
Se sappiamo che una certa successione an è geometrica, allora la conoscenza di
due qualsiasi dei suoi termini permette di determinare l’intera successione. Per ogni
coppia di numeri naturali distinti h, k si ha infatti:
ah
a0 q h
=
= q h−k
ak
a0 q k
e quindi
r
q=
h−k
ah
ak
e
a0 =
ah
.
qh
S. Console – M. Roggero
Le successioni
119
Esempio
1.11. La successione geometrica an tale che a3 = 8 e a7 = 32 ha ragione
q
√
√
√
√ √
a7
3
4
q = a3 = 4 4 = 2 e quindi a0 = (√a2)
2 ossia è data da an = 2 2·( 2)n =
3 = 2
√
( 2)n+3 .
La formula che dà la somma dei primi n termini di una successione geometrica
(di ragione q 6= 1) è:
1 − q n+1
a0 + a1 + · · · + an = a0
1−q
§ A. 1.3
Applicazioni delle successioni
Le successioni aritmetiche e geometriche intervengono nello studio di numerosi problemi di tipo economico, biologico, medico.
Esempio 1.12. Si vuole trovare una formula che dia il valore dello stipendio di un
lavoratore dopo n anni, sapendone il valore iniziale S0 e supponendone un aumento
annuale pari al 2% di S0 . Procedendo ricorsivamente, abbiamo:
stip. iniziale: S(0) = S0
stip. dopo 1 anno: S(1) = S0 +
2
100 S0
stip. dopo 2 anni: S(2) = S(1) +
......
2
100 S0
stip. dopo n anni: S(n) = S(n − 1) +
= (S0 +
2
100 S0
2
100 S0 )
+
2
100 S0
2
= S0 + 2( 100
S0 )
2
= S0 + n 100
S0 .
Il problema è descritto da una successione aritmetica di termine iniziale S0 e ragione
2
100 S0 : in tal caso si dice anche che lo stipendio cresce in modo lineare.
Esempio 1.13. Si vuole trovare una formula che dia il valore dello stipendio Q(n)
di un lavoratore dopo n anni, sapendone il valore iniziale S0 e supponendone un
aumento annuale pari al 2% dello stipendio dell’anno precedente. Procedendo ricorsivamente, abbiamo:
stip. iniziale: Q(0) = S0
stip. dopo 1 anno: Q(1) = S0 +
2
100 S0
stip. dopo 2 anni: Q(2) = Q(1) +
......
=
2
100 Q1
102
100 S0
=
102
100 Q1
=
102 2
S0
100
Modelli Matematici applicati all’Ecologia (03.04.06)
120
Appendice 1
stip. dopo n anni: Q(n) = Q(n − 1) +
2
100 Qn−1
=
102
100 Qn−1
=
102 n
S0 .
100
Il problema è descritto dà una successione geometrica di termine iniziale S0 e ragione
q = 102
100 : in tal caso si dice anche che lo stipendio cresce in modo esponenziale.
Esempio 1.14. Si vuole schematizzare in modo ricorsivo il processo di decadimento radioattivo. Alcune sostanze decadono nel tempo, trasformandosi in altre
sostanze; si dice tempo di dimezzamento il periodo T in cui decade la metà degli
atomi. Assumendo come unità di misura dei tempi T e indicando con N0 il numero
degli atomi presenti inizialmente si ha:
numero iniziale: N (0) = N0
dopo 1 periodo: N (1) = 21 N0
dopo 2 periodi: N (2) = 12 (N (1)) = 14 N0
......
dopo n periodi: N (n) = 12 N (n − 1) = ( 12 )n N0 .
Il processo è descritto dà una successione geometrica di termine iniziale N e ragione
1
2.
Esempio 1.15. In condizioni ideali la crescita di una popolazione di batteri
ha un andamento di tipo esponenziale. In altre parole, se rileviamo il numero di
batteri presenti in quella popolazione e ripetiamo il conteggio a distanza regolare
di tempo, i numeri ottenuti formano una successione geometrica an = Bq n dove
a0 = B è il numero di batteri inizialmente presenti e q dipende dal tipo di batteri e
dalle condizioni ambientali.
Il tempo T necessario perché il numero di batteri di quella popolazione raddoppi
si dice appunto tempo di raddoppio. Assumendo come intervallo di tempo tra una
rilevazione e l’altra proprio T , la successione assume la forma
an = 2n B
dove an è il numero di batteri presenti dopo un periodo nT pari a n tempi di
raddoppio.
Più in generale, se l’intervallo tra due rilevazioni è un qualsiasi tempo t, allora
la successione an che fornisce il numero di batteri presenti dopo n rilevazioni, ossia
dopo un tempo nt è
t
an = q n B dove la ragione è q = 2( T ) .
In modo più semplice, se t ∼ νT (oppure νT < t < (ν + 1)T ) per un opportuno
numero intero ν, allora dopo un tempo t il numero Bt dei batteri presenti nella
popolazione sarà Bt ∼ 2ν B (rispettivamente 2ν B < Bt < 2ν+1 B).
S. Console – M. Roggero
Le successioni
121
Esempio 1.16. Un modello piu’ aderente alla realtà per la crescita di una popolazione è dato da successioni del tipo xn+1 = f (xn ) con f (x) = x(2 − x) (modello
logistico).
Per visualizzare i termini della successione
si possono disegnare la
funzione f (x) e la bisettrice del 1 e 3 quadrante.
Usando la diagonale si
possono allora visualizzare i succesivi termini
(basta proiettare le linee verticali sull’asse x)
§ A. 1.4
1
1
0.8
0.8
0.6
0.6
y
0.4
0.4
0.2
0
0.2
0.2
0.4
t
0.6
0.8
1
0
0.2
0.4
t
0.6
0.8
1
Esercizi risolti
A.1 Il tempo di dimezzamento del C14 è di circa 5730 anni.
(a) Determinare l’età approssimativa di un reperto fossile nel quale la concentrazione di C14 risulta
il 12% di quella dell’analogo organismo vivente.
(b) Determinare la concentrazione di C14 in un reperto fossile di circa 23˙000 anni.
12
Soluzione: (a) Poiché 12% = 100
= 0, 12 ∼ ( 12 )3 , il reperto in questione avrà un’età pari a circa 3
periodi di dimezzamento, ossia 5730 · 3 = 17.190 anni.
(b) Poiché 23˙000 ∼ 5730 · 4, la concentrazione nel reperto fossile sarà pari a ( 21 )4 =
1
16
∼ 7% della
concentrazione iniziale.
A.2 In una coltura batterica ci sono inizialmente N batteri, che raddoppiano ogni 1600 .
(a) Quanti batteri ci saranno dopo 22 h?
(b) Dopo quanti minuti ci sono nella coltura il 25% del numero finale di batteri trovato al punto
(a)?
Soluzione: (a) 22h = 13200 , ossia circa 8 volte il tempo di raddoppio; quindi dopo 22 ore troveremo
circa 28 N = 256N batteri.
(b) Il 25% di 256 è 64 = 26 e 6 tempi di raddoppio sono 960 minuti.
A.3
Il tempo di dimezzamento dell’ossigeno O15 è 12400 .
(a) Indicata con Q la concentrazione iniziale, trovare quella che si ha dopo 90 .
(b) Dire dopo quanti secondi si ha il 25% di Q.
Soluzione: (a) Poiché 90 = 54000 e si ha 124 · 4 < 540 < 124 · 5, allora in 90 l’O15 si dimezza tra 4
e 5 volte; quindi la concentrazione sarà compresa tra ( 12 )4 Q e ( 12 )5 Q.
(b) Poiché
25
100
= ( 21 )2 , si ha il 25% di Q dopo due tempi di dimezzamento, quindi dopo circa 24800 .
A.4 Ho una sostanza con delle impurità ed uso un procedimento che ad ogni applicazione elimina
il 30% delle impurità. Se uso il procedimento 5 volte, qual è la percentuale di impurità eliminate?
Modelli Matematici applicati all’Ecologia (03.04.06)
122
Appendice 1
7
Soluzione: Dopo ogni procedimento resta il 10
= 70% delle impurità che c’erano prima. Perciò
5
0.168
=
16.8%
delle impurità che c’erano all’inizio. Sono
dopo 5 procedimenti resta il ( 7 ) ∼
=
10
eliminate quindi circa 1 − 0.168 = 0.832 = 83.2% delle impurità.
A.5 Tenendo conto che il tempo di dimezzamento del 14 C è di 5730 anni, determinare con ragionevole approssimazione quale età ha un reperto fossile la cui concentrazione di 14 C è il 3% di quella
originaria.
Soluzione: La concentrazione Cn di 14 C dopo n dimezzamenti è
1
Cn = C0 ( )n ,
2
dove C0 è la concentrazione iniziale. Dunque se la concentrazione dopo n dimezzamenti è il 3% di
quella originaria, si ha:
Cn
3
1
=
∼ ( )5 .
C0
100
2
Quindi n = 5 e gli anni sono 5 × 5730 = 28650.
A.6 Calcolare il 3% del 10% di una quantità a.
3
10
Soluzione: Moltiplicando a per 100
e per 100
otteniamo
30
a
10000
=
0,3
100
ossia lo 0, 3% di a.
Si faccia attenzione a non calcolare il prodotto di 3% per 10% come (3 · 10)% = 30%!!
A.7 La rendita catastale di un alloggio, rivalutata del 5%, vale 525e. Calcolare la rendita effettiva
(senza tale rivalutazione).
Soluzione: indichiamo con R il valore cercato. Sappiamo che R + 5%R = 525. Risolvendo
l’equazione
105
R
100
= 525, troviamo R = 500e.
A.8 Ricordiamo che per concentrazione di una soluzione si intende il rapporto tra la quantità di
soluto e la quantità di soluzione. Per esempio, dire che una soluzione ha concentrazione pari al
7% significa che vi sono 7 grammi di soluto (e quindi 93 grammi di solvente) ogni 100 grammi di
soluzione. Supponiamo ora di avere 5 chilogrammi di una soluzione concentrata al 40%. Si chiede
quanto solvente occorra aggiungere per ottenere una soluzione concentrata al 25%.
Soluzione: 5 Kg di soluzione al 40% sono costituiti da 2 Kg di soluto e 3 Kg di solvente. Vogliamo
una nuova soluzione S di cui i 2 Kg di soluto siano il 25% del totale; allora S è tale che
2
S
=
25
100
= 14 .
Otteniamo S = 8 Kg; vi saranno quindi 6 = 8 − 2 Kg di solvente anziché i 3 kg iniziali. Occorre
quindi aggiungere 3 Kg di solvente per diluire al 25% la concentrazione della soluzione.
A.9 All’esame di Matematica si presenta il 70% degli studenti che si erano iscritti. Durante lo
svolgimento della prova scritta si ritira il 40% degli studenti presenti. Supera lo scritto il 50% degli
studenti che lo hanno consegnato. Trovare la percentuale di studenti iscritti che hanno superato la
prova scritta.
Soluzione: indichiamo con I il numero degli studenti iscritti all’esame e riscriviamo i dati:
70
I
100
60
·
100
50
·
100
è il numero di studenti presenti alla prova;
70
I
100
42
I
100
=
=
42
I
100
21
I
100
è il numero degli studenti non ritirati e quindi dei compiti consegnati;
è il numero di studenti che superano lo scritto.
Superano quindi l’esame il 21% degli studenti iscritti.
A.10 Una azienda ha prodotto nel 2004 il 50% in più rispetto al 2003 e nel 2005 il 50% in meno
del 2004. Esprimere in percentuale quanto la stessa azienda ha prodotto nel 2005 rispetto al 2003.
Soluzione: indichiamo con P3 , P4 , P5 , i prodotti degli anni 2003, 2004 e 2005 rispettivamente.
S. Console – M. Roggero
Le successioni
50
P = 150
P .
100 3
100 3
50
50 150
50
P
=
P
= 100
( 100 P3 )
100 4
100 4
123
Abbiamo P4 = P3 +
Quindi P5 = P4 −
=
75
P .
100 3
A.11 Definiamo prevalenza di una malattia il rapporto, a una dato istante, tra i soggetti malati e
l’intera popolazione in esame. Supponiamo che gli abitanti dell’Italia siano circa 57.400.000, quelli
del nord circa 25.500.000, quelli delle isole 7.000.000 circa. Sapendo che una certa malattia ha una
prevalenza dello 0, 7% al nord, è assente al centro e al sud, ha una prevalenza dell’1, 3% nelle isole,
calcolare in percentuale la prevalenza della malattia sulla popolazione italiana.
Soluzione: lo 0, 7% di 25.500.000 vale 178.500, l’1, 3% di 7.000.000 è 91.000. Il rapporto tra la
somma di questi numeri, pari a 269.500, e 57.400.000 è circa 0,0047. La percentuale richiesta è
dunque 0, 47% circa.
§ A. 1.5
Altri esercizi
A.12 La costruzione di una diga su un torrente ha dato luogo alla formazione di un laghetto della
capienza di circa 12˙000m3 d’acqua. Ma il torrente sbarrato immette mediamente nel bacino 200
m3 di detriti l’anno. Determinare la capienza del laghetto nei primi 20 anni mediante una formula.
Dopo quanti anni la capienza del laghetto sarà ridotta alla meta di quella iniziale?
A.13 Determinare il termine generale Sn di una successione geometrica tale che S0 = 20 e S10 = 60;
A.14 Provare che la successione Cn =
√
n non è nè aritmetica , nè geometrica.
A.15 Determinare ragione e valore iniziale della successione geometrica Sn tale che S5 = 21 e
S6 = 14.
A.16 Stabilire se i numeri 81, 27, 9, 3 sono i termini iniziali di una progressione aritmetica e/o
geometrica e, in caso affermativo, determinare i due termini successivi.
A.17 Una colonia di batteri raddoppia ogni 10 giorni. Se dopo 10 settimane abbiamo 20 milioni di
batteri, quanti ne avevamo all’inizio?
A.18 Nella fase di lavaggio di una lavatrice sono immessi 90g di detersivo. In ogni risciacquo viene
eliminato il 96% di detersivo presente in quel momento. Quanti risciacqui sono necessari affinché
rimanga meno di 0.1 g di detersivo?
A.19 Supponiamo che la legge che regola la crescita nel tempo di una colonia di microorganismi
(controllati con una cadenza di un giorno) sia una progressione geometrica. Se dopo 5 giorni i
microorganismi sono 9˙000 e dopo 8 giorni 72˙000, qual è la ragione della progressione geometrica?
Dopo quanto tempo i microorganismi sono 1 milione?
A.20 Una colonia di microorganismi raddoppia ogni 14 giorni. Se dopo 40 settimane ne abbiamo
30 milioni, quanti microorganismi avevamo all’inizio?
A.21 Si consideri una successione che ha i seguenti due termini: A2 = 3, A7 = 15.
Modelli Matematici applicati all’Ecologia (03.04.06)
124
Appendice 1
a) Scrivere il termine generale nel caso in cui si tratti di una successione aritmetica.
b) Scrivere il termine generale nel caso in cui si tratti di una successione geometrica.
A.22 In un vivaio ogni anno muoiono circa il 25% delle piante presenti. Da quante nuove piantine
bisogna partire per avere 300 piante di 5 anni?
A.23 (a) Determinare l’età approssimativa di un reperto fossile nel quale la concentrazione di C14
risulta il 26% di quella dell’analogo organismo vivente .
(b) Determinare la concentrazione di C14 in un reperto fossile di circa 29˙000 anni.
S. Console – M. Roggero
Appendice 2
I numeri complessi
§ A. 2.1
Definizione di C
Introduciamo un numero nuovo, che chiameremo i che sia una soluzione dell’equazione X 2 + 1 = 0, ossia che abbia la proprietà che il suo quadrato fa −1:
i2 = −1.
Affinchè questa costruzione abbia un senso, non possiamo limitarci ad ampliare
i numeri reali mediante la sola aggiunta del nuovo numero i, ma dobbiamo far sı̀ che
si possano eseguire anche somme e prodotti.
Dovremo quindi introdurre anche tutti
√
2
i numeri del tipo 2i, 3 + i, (i + 5)( 2i − π) e cosı̀ via.
Definizione 2.1. Chiamiamo numeri complessi tutte le espressioni della forma
a + ib dove a e b sono numeri reali.
Il numero complesso i si dice unità immaginaria.
Se z = a+ib, a si dice parte reale di z, denotata Re(z), e b si dice coefficiente
dell’immaginario di z, denotato Im(z).
Due numeri complessi sono uguali se e solo se coincidono sia la parte
reale sia quella immaginaria.
Definiamo quindi le operazioni di somma e prodotto in C a partire dalle operazioni di R secondo le regole del calcolo algebrico nel modo seguente:
(a + ib) + (c + id) =
a + ib + c + id
= (a + c) + i(b + d)
(a + ib) · (c + id) = ac + iad + ibc + i2 bd = ac + iad + ibc − bd
= (ac − bd) + i(ad + bc)
Per come sono state definite le operazioni, possiamo pensare ogni numero reale
a come un particolare numero complesso che ha parte immaginaria nulla, ossia a è
anche il numero complesso a + i0.
In C valgono tante delle proprietà delle operazioni che conosciamo in R. Eccone
alcune:
125
126
Appendice 2
1) la somma e il prodotto sono associative e commutative e vale la proprietà distributiva del prodotto rispetto alla somma;
2) c’è l’elemento neutro rispetto alla somma 0 = 0 + i0, infatti 0 + (a + ib) = a + ib;
3) c’è l’opposto di ogni elemento −(a + ib) = (−a) + i(−b);
4) c’è l’elemento neutro rispetto al prodotto 1 = 1 + i0, infatti 1 · (a + ib) = a + ib;
5) ogni numero complesso, tranne 0, ha un inverso (a + ib)−1 =
a
a2 +b2
b
− i a2 +b
2.
Dalla validità di 1), 2) 3), 4), 5) segue che C è un campo.
Ritorniamo ora sull’equazione x2 + 1 = 0 per osservare che in C risulta avere,
oltre alla soluzione x = i, anche la soluzione x = −i, come si può verificare
immediatamente tramite sostituzione.
Il numero complesso −i ha in realtà le stesse proprietà di i e risulta essere
perfettamente interscambiabile con i.
Se z = a + ib è un qualsiasi numero complesso, scambiando in esso i con −i
otteniamo il numero complesso z = a − ib, che si dice coniugato di z.
Valgono le seguenti proprietà relative al coniugato di un numero complesso:
- la somma di un numero e del suo coniugato: z + z è sempre un numero reale;
- i numeri reali sono gli unici numeri complessi che coincidono col loro coniugato
z = z ⇐⇒ z ∈ R;
- il prodotto di un numero diverso da 0 e del suo coniugato z · z è sempre un numero
reale strettamente positivo.
- il coniugato di una somma coincide con la somma dei coniugati, ossia
z1 + z 2 = z 1 + z2
- il coniugato di un prodotto coincide col prodotto dei coniugati, ossia
z1 · z 2 = z 1 · z2
Si dice modulo del numero complesso z = a + ib il numero reale positivo o nullo
p
√
|z| = z · z = a2 + b2 .
Da quanto sopra risulta |z| = 0 se e solo se z = 0.
S. Console – M. Roggero
I numeri complessi
127
Possiamo allora rivedere in altro modo la formula dell’inverso di un numero
complesso:
z −1 =
§ A. 2.2
z
.
z·z
Il teorema fondamentale dell’algebra
In C tutte le equazioni di grado 2 a coefficienti reali x2 + px + q = 0 hanno due
soluzioni complesse, eventualmente coincidenti.
− p2
Se √infatti il discriminante
∆ = p2 −4q è positivo o nullo, ci sono le due radici reali
√
p
∆
∆
+ 2 e − 2 − 2 , mentre se ∆ è negativo ( e quindi −∆ è un numero reale
√
√
e − p2 − i −∆
positivo!), ci sono le due radici complesse (non reali) − p2 + i −∆
2
2 .
Si noti che queste ultime sono due numeri complessi coniugati ossia ottenibili
uno dall’altro scambiando i con −i ossia passando ai coniugati!
Introducendo semplicemente il nuovo numero i siamo cosı̀ in grado di risolvere
tutte le equazioni di secondo grado. Non facile da provare, ma ugualmente vero è il
seguente, più generale, risultato.
Teorema fondamentale dell’algebra: Ogni equazione polinomiale di grado
positivo n a coefficienti reali o complessi
xn + a1 xn−l + a2 xn−2 + ... + an−1 x + an = 0
ha sempre n soluzioni complesse (eventualmente coincidenti) α1 , α2 , . . . , αn , e quindi
il polinomio
F (x) = xn + a1 xn−l + a2 xn−2 + ... + an−1 x + an
si decompone nel prodotto di polinomi di grado 1:
F (x) = (x − α1 )(x − α2 ) . . . (x − αn ).
Se l’equazione ha coefficienti reali e se α è una sua soluzione, allora anche il
coniugato α è soluzione e il polinomio si decompone nel prodotto di fattori di grado
1 corrsipondenti alle soluzioni reali e di fattori di grado due corrispondenti alle coppie
di soluzioni complesse coniugate.
§ A. 2.3
Esercizi
A.1 Ridurre le espressioni alla forma a + ib specificando qual è la parte reale e quale immaginaria
Modelli Matematici applicati all’Ecologia (09.12.06)
128
Appendice 2
1. (3 + 4i) + (2 + 5i);
2. (3 + 4i)(4 + 5i);
3. (1 + 3i)(6 − 6i);
4. (2 + 3i)(2 − 3i);
5. 3i(4 + i)2
A.2 Calcolare tutte le soluzioni delle equazioni:
√
√
√
√
1. (x − 3)(4 − 5x)( 3 − 5x)( 3 + 5x) = 0
√ 4
2. (x − 2) = 0
3. (x2 − 3x − 5)(3x − 2x2 + 1) = 0
4. x2 − 5x + 7 = 0
√
√
A.3 Verificare che i numeri complessi (1 + i 3) e (1 − i 3) (oltre che naturalmente 2) sono radici
del polinomio x3 + 8. Eseguire la verifica sia mediante sostituzione, sia applicando il teorema di
Ruffini.
A.4 Calcolare la parte reale e la parte immaginaria dei seguenti numeri complessi:
(1 + i)(2i − 3),
(2 − i)2 ,
(1 − 3i)3 ,
6 + 5i
,
3−i
3 − 2i
2 − 3i
+
.
1 + 5i
2−i
A.5 Dire se il numero complesso −1+i è oppure no una soluzione dell’equazione x3 −2x2 +5x−2 = 0.
A.6 Verificare che il coniugato di una somma è la somma dei coniugati e che il coniugato di un
prodotto è il prodotto dei coniugati.
A.7 Sia z un qualsiasi numero complesso. Verificare che z + z e z · z sono numeri reali e che z − z
è immaginario puro.
A.8 Sia z un numero complesso.Verificare che z = z se e solo se z ∈ R.
A.9 Cosideriamo il polinomio F (X) = 2X 5 − 13X 4 + 37X 3 − 57X 2 + 48X − 18.
Verificare che 1 − i è radice di F (X).
Trovare tutte le radici razionali di F (X).
Determinare la fattorizzazione di F (X) in fattori irriducibili in R[X].
A.10 Calcolare a mano (se possibile) le soluzioni delle seguenti equazioni e poi risolverle mediante
Derive:
a. x3 − 1 = 0
b. x4 − 1 = 0
c. x4 + 1 = 0
d. x3 − 3x + 1 = 0
5
e. x − 3x + 1 = 0.
S. Console – M. Roggero
Appendice 3
Le funzioni di due o più variabili
§ A. 3.1
grafica
Funzioni a più variabili e loro rappresentazione
Indichiamo con R2 , l’insieme delle coppie di numeri reali, cioè i punti del piano
euclideo P = (x, y). Una funzione di 2 variabili è una legge (una regola, una
formula) che associa un numero z ad ogni punto P ∈ R2
z = f (x, y) = f (P ).
Una tale funzione si dice, più precisamente, funzione reale (o a valori reali) nelle
2 variabili reali (x, y). Queste ultime si dicono variabili indipendenti, mentre
la z si dice variabile dipendente. In generale, la funzione f non è definita in tutti
i punti P di R2 , ma su un sottoinsieme D (dominio).
Si può generalizzare questo concetto a quello di funzioni di n variabili. (Nel caso
n = 1 ritroviamo il concetto di funzione ad una variabile y = f (x).)
Per le funzioni a due variabili si può
dare una rappresentazione grafica, detta anche rappresentazione topografica,
nello spazio tridimensionale euclideo, riferito ad una terna di assi cartesiani (x, y, z)
con origine in un punto O. Il grafico di
una funzione z = f (x, y) è costituito da
quei punti P = (x, y, z) aventi come prime
due coordinate i valori di (x, y) ∈ D e come
terza coordinata il valore z = f (x, y). In
genere (cioè per funzioni non troppo complicate) il grafico è una superficie, detta superficie topografica, nello spazio tridimensionale che si proietta sopra il dominio
D del piano (x, y). Questa superficie ha la proprietà di incontrare in un solo punto
le parallele all’asse z uscenti dai punti di D.
129
130
§ A. 3.2
Appendice 3
Derivate parziali
Data una funzione z = f (x, y), si può pensare di tener fissa la variabile y (considerandola quindi come una costante) e di derivare la funzione f (x, y) come funzione della
sola variabile x. Si ottiene il questo modo (ammesso che esista) la derivata parziale
di f (x, y) rispetto a x. La si denota in vari modi:
∂f
,
∂x
fx0 .
Analoga operazione si può fare considerando x costante; si ottiene la derivata parziale
rispetto a y, denotata con
∂f
,
fy0 .
∂y
Per esempio, le due derivate parziali della funzione z = f (x) = sin(3x + y 2 ) sono
rispettivamente
∂f
= fx0 = 3 cos(3x + y 2 ),
∂x
∂f
= fy0 = 2 y cos(3x + y 2 ).
∂y
Va osservato che le derivate parziali sono ancora funzioni delle due variabili (x, y),
per cui ha senso procedere ad ulteriori derivazioni. Si giunge cosı̀ (ammesso che
esistano) alle derivate parziali seconde che si denotano con
∂2f
,
∂x2
∂2f
,
∂y∂x
∂2f
,
∂x∂y
∂2f
,
∂y 2
oppure con
00
fxx
,
00
fxy
,
00
fyx
,
00
fyy
,
Nell’esempio considerato abbiamo
00
fxx
=
− 9 sin(3x + y 2 ),
00
fxy
=
− 6 y sin(3x + y 2 ),
00
fyx
=
− 6 y sin(3x + y 2 ),
00
fyy
= 2 cos(3x + y 2 ) − 4 y 2 sin(3x + y 2 ).
00 e f 00 coinSi osserva da quest’esempio che le due derivate parziali miste fxy
yx
cidono. Non si tratta di un caso particolare ma di una notevole proprietà generale,
valida sotto certe ipotesi sulla funzione e sulle sue derivate prime e seconde, di solito
soddisfatte nelle normali applicazioni (per esempio l’esistenza e la continuità delle
derivate seconde). Questa proprietà va sotto il nome di teorema d’inversione
dell’ordine delle derivazioni o di teorema di Schwarz:
∂2f
∂2f
=
.
∂y∂x
∂x∂y
S. Console – M. Roggero
Le funzioni di due o più variabili
131
L’operazione di derivata parziale si può reiterare ottenendo, se esistono, derivate
parziali terze, quarte, ecc. Vale ancora per le derivate miste e subordinatamente alle
ipotesi accennate, il teorema d’inversione dell’ordine delle derivazioni, il quale allora
va interpretato con il fatto che gli operatori di derivata parziale
∂
,
∂x
∂
,
∂y
commutano, possono cioè applicarsi in ordine qualsiasi.
Modelli Matematici applicati all’Ecologia (09.12.06)
Bibliografia
[BC]
F. Brauer, C. Castillo-Chavez: Mathematical Models in Population Biology
and Epidemiology, Springer , TAM 40.
[C]
A. C. Capelo, Modelli matematici in biologia, introduzione all’ecologia
matematica, Zanichelli Decibel editore.
[F]
J. C. Frauenthal: Introduction to Population Modeling, Birkhäuser,
Boston, Basel, Stuttgart, 1980.
[ML]
Marsili-Libelli: Modelli Matematici in Ecologia, Pitagora
132