6.1 Συνεχή Μοντέλα Αλληλεπιδρώντων Πληθυσμών

Το πιο απλό μοντέλο θηράματος–θηρευτή (prey – predator), αλλά και συγχρόνως ίσως το πιο «διάσημο» μοντέλο Μαθηματικής Βιολογίας, είναι το μοντέλο Lotka-Volterra, το οποίο αναπτύχθηκε από τον Volterra το 1926 [Volterra, V., “Variazioni e fluttuazioni del numero d’individui in specie animali conviventi”, Mem. Acad. Lincei Roma, 2, 31–113, (1926), για να εξηγήσει τους περιοδικά ταλαντωμένους πληθυσμούς δύο ειδών ψαριών στην Αδριατική. Βασικές υποθέσεις του μοντέλου αυτού είναι ότι: (α) ο θηρευτής δεν μπορεί να τραφεί παρά μόνο με το θήραμα, ενώ (β) για το θήραμα υπάρχει «απεριόριστη» τροφή στο περιβάλλον. Αρχικά θα μελετήσουμε τη συνεχή και κατόπιν τη διακριτή του μορφή. [26, 37]

Έστω δύο είδη ψαριών από τα οποία το ένα τρέφεται με το άλλο και, οι πληθυσμοί των δύο ειδών περιγράφονται από συνεχείς μεταβλητές ως προς το χρόνο

• •

  $N(t)$: ο πληθυσμός του θηράματος (prey), και

• •

  $P(t)$: ο πληθυσμός του θηρευτή (predator), σε χρόνο $t$.

Στόχος μας είναι η δημιουργία κατάλληλων μαθηματικών σχέσεων (στην περίπτωσή μας συνήθων διαφορικών εξισώσεων), οι οποίες περιγράφουν τη χρονική εξέλιξη των δύο πληθυσμών. Ξεκινάμε με συγκεκριμένα (βιολογικές) υποθέσεις. Θεωρούμε ότι:

1. 1.

   Ο πληθυσμός του θηράματος, απουσία του θηρευτή μεγαλώνει εκθετικά. Αυτό σημαίνει ότι στο ρυθμό μεταβολής (παράγωγος) του πληθυσμού του θηράματος θα υπάρχει θετικός όρος ($+aN$).

2. 2.

   Αντίστοιχα ο πληθυσμός του θηρευτή, απουσία του θηρευτή μειώνεται εκθετικά. Αυτό σημαίνει ότι στο ρυθμό μεταβολής (παράγωγος) του πληθυσμού του θηρευτή θα υπάρχει αρνητικός όρος ($-dN$).

3. 3.

   Η μείωση του πληθυσμού του θηράματος λόγω του θηρευτή είναι ανάλογη των «επαφών» θηράματος-θηρευτή, δηλαδή και μάλιστα είναι γραμμικής μορφής ως προς το γινόμενο $NP$. Δηλαδή στο ρυθμό μεταβολής του πληθυσμού του θηράματος υπάρχει αρνητικός όρος ($-bNP$).

4. 4.

   Αντίστοιχα η αύξηση του πληθυσμού του θηρευτή λόγω του θηράματος είναι ανάλογη των «επαφών» θηράματος-θηρευτή, δηλαδή και μάλιστα είναι γραμμικής μορφής ως προς το γινόμενο $NP$. Δηλαδή στο ρυθμό μεταβολής του πληθυσμού του θηράματος υπάρχει θετικός όρος ($+cNP$).

Στους παραπάνω όρους οι $a,b,c,d$ είναι θετικές σταθερές. Τελικά το συνεχές μοντέλο Lotka-Volterra περιγράφετε ως ένα σύστημα δύο (συνήθων) διαφορικών εξισώσεων:

$$\begin{array}{c}\hfill \frac{dN}{dt}=N(a-bP)\\ \hfill \frac{dP}{dt}=P(cN-d).\end{array}$$

(6.1)

Θα πρέπει εδώ να σημειώσουμε ότι αν και το παραπάνω μοντέλο (σύστημα εξισώσεων) εμφανίστηκε όπως αναφέραμε για το συγκεκριμένο βιολογικό πρόβλημα το 1926, είχε αναπτυχθεί ανεξάρτητα από τον Lotka το 1920 για να περιγράψει μια χημική αντίδραση η οποία εμφανίζει περιοδική συμπεριφορά στις συγκεντρώσεις των στοιχείων που συμμετέχουν σε αυτή! Για αυτό το λόγο ονομάζεται Lotka-Volterra μοντέλο. Εύκολα καταλαβαίνουμε ότι το κοινό στοιχείο των παραπάνω προβλημάτων (η περιοδικότητα στους πληθυσμούς ή στις συγκεντρώσεις των στοιχείων) θα πρέπει να εμφανίζεται λόγω της πεπλεγμένης δομής των συγκεκριμένων εξισώσεων.

Ξεκινάμε τη μελέτη του παραπάνω μοντέλου κάνοντας μια κατάλληλη αδιαστατοποίηση ώστε να «απλοποιήσουμε» τη μορφή των εξισώσεων. Η αδιαστατοποίηση των μεταβλητών είναι συχνά μια πολύ χρήσιμη μαθηματική διαδικασία όπως είδαμε αναλυτικά στο Κεφάλαιο 1. Για τις παραπάνω εξισώσεις επιλέγουμε τις εξής αδιάστατες μεταβλητές.

$$u(\tau )=\frac{cN(t)}{d},v(\tau )=\frac{bP(t)}{a},\tau =at,\tilde{d}=\frac{d}{a}.$$

(6.2)

Με την επιλογή αυτή το αρχικό σύστημα εξισώσεων γίνεται (η απόδειξη αφήνεται ως άσκηση):

(6.3)

Το παραπάνω σύστημα εξισώσεων αποτελεί την αδιάστατη μορφή του (συνεχούς) μοντέλου Lotka-Volterra. Τα σημεία ισορροπίας του συστήματος, είναι οι τιμές των μεταβλητών $u,v$ για τις οποίες

$${\frac{du}{d\tau }|}_{u^{*},v^{*}}={\frac{dv}{d\tau }|}_{u^{*},v^{*}}=0$$

Το παραπάνω σύστημα έχει δύο σημεία ισορροπίας (fixed points):

• Α)

  $(u^{*},v^{*})=(0,0)$

• Β)

  $(u^{*},v^{*})=(1,1)$

Κατόπιν εξετάζουμε την ευστάθεια του μοντέλου, ορίζοντας μικρή διαταραχή $(x,y)$ στον πληθυσμό $(u,v)$ και κάνοντας την αντίστοιχη γραμμικοποίηση, δηλαδή ανάπτυγμα Taylor και κράτημα μόνο των γραμμικών όρων.

Α) Για το πρώτο σημείο $(u^{*},v^{*})=(0,0)$ έχουμε τη διαταραχή:

	$u$	$=u^{*}+x=x$
	$v$	$=v^{*}+y=y$

Με γραμμικοποίηση θα καταλήξουμε (η απόδειξη αφήνεται ως άσκηση στον αναγνώστη), για τη χρονική εξέλιξη της διαταραχής, στις εξισώσεις:

$$\frac{dx}{d\tau }\approx x,\frac{dy}{d\tau }\approx -\tilde{d}y.$$

(6.4)

Λύσεις του παραπάνω συστήματος είναι:

$$\left[\begin{array}{c}\hfill x(\tau )\hfill \\ \hfill y(\tau )\hfill \end{array}\right]=B{e}^{\lambda \tau }$$

(6.5)

με ιδιοτιμές:

Συνεπώς η ευστάθεια καθορίζεται ως:

• •

  Αν ${\lambda }_1>0$ τότε η διαταραχές $x(\tau )$ και $y(\tau )$ μεγαλώνουν εκθετικά, δηλαδή το σημείο $(u,v)=(0,0)$ είναι ασταθές.

• •

  Αν ${\lambda }_1>0$ και τότε και το σημείο είναι σαγματικό (saddle point).

Β) Για το δεύτερο σημείο $(u^{*},v^{*})=(1,1)$ ορίζουμε μικρή διαταραχή

	$u$	$=u^{*}+x=1+x$
	$v$	$=v^{*}+y=1+y$

Αντίστοιχη γραμμικοποίηση δίνει:

	${\displaystyle \frac{d(x+1)}{d\tau }}$	$\approx (1+x)(1-1-y)=-y(1+x)=-y-xy\approx -y$
	${\displaystyle \frac{d(y+1)}{d\tau }}$	$\approx \tilde{d}(1+y)(1+x-1)=\tilde{d}x+\tilde{d}xy\approx \tilde{d}x$

εφόσον ισχύει: $xy\ll y$ και $xy\ll x$.

Οι ιδιοτιμές του παραπάνω συστήματος είναι:

δηλαδή το σημείο $(u,v)=(1,1)$ στο γραμμικοποιημένο σύστημα είναι κέντρο. Η παραπάνω λύση δίνει τα παρακάτω δύο ιδιοδιανύσματα, $l,m$:

$$\left[\begin{array}{c}\hfill x(\tau )\hfill \\ \hfill y(\tau )\hfill \end{array}\right]=l\exp \mathrm{}\left(i\sqrt{\tilde{d}}\tau \right)+m\exp \mathrm{}\left(-i\sqrt{\tilde{d}}\tau \right)$$

(6.6)

Βασικό πρόβλημα του μοντέλου Lotka-Volterra είναι ότι οι λύσεις δεν είναι σταθερές (έλλειψη) με τη συνήθη έννοια, δηλαδή όπως είδαμε μια μικρή διαταραχή οδηγεί σε άλλη λύση. Αυτό είναι χαρακτηριστικό των συστημάτων διατήρησης. Παρόλα αυτά έχει εφαρμοστεί με σχετική επιτυχία σε διαφορετικά πραγματικά συστήματα.

Η λύση του αρχικού συστήματος μπορεί να γίνει με αριθμητικές μεθόδους. Μια τέτοια λύση πραγματοποιήθηκε, σε γλώσσα Matlab, χρησιμοποιώντας μεθόδους Runge-Kutta (ο κώδικας υπάρχει στο τέλος του κεφαλαίου). Στα Σχήματα 6.1α και 6.1β παρουσιάζεται η εξέλιξη των πληθυσμών του θηρευτή και του θηράματος για δυο διαφορετικούς συνδυασμούς τιμές των παραμέτρων του προβλήματος, δηλαδή των σταθερών $a,b,c,d$ και για αρχικές τιμές των δύο πληθυσμών $N(t=0)=P(t=0)=1.0$.

Όπως βλέπουμε η συμπεριφορά και των δύο πληθυσμών είναι περιοδική στο χρόνο.

Πιο συγκεκριμένα στο Σχήμα 6.1α βλέπουμε την συμπεριφορά των πληθυσμών για τιμές των παραμέτρων a=1.0, b=2.0, c=2.0, d=1.0. Στην περίπτωση αυτή οι δύο σταθερές $b$ και $c$ που αφορούν την μείωση του πληθυσμού του θηράματος, και την αντίστοιχη αύξηση του πληθυσμού του θηρευτή, λόγω της αλληλεπίδραση θηράματος – θηρευτή, είναι ίσες. Επίσης και οι δύο σταθερές που περιγράφουν την αύξηση (μείωση) του πληθυσμού του θηράματος (θηρευτή) απουσία του θηρευτή (θηράματος) είναι ίσες. Για αυτό το λόγο η συμπεριφορά των δύο πληθυσμών είναι «παρόμοια», δηλαδή η «περίοδος» και το «πλάτος» της ταλάντωσης είναι το ίδιο.

Στο Σχήμα 6.1β βλέπουμε την συμπεριφορά των πληθυσμών για τιμές των παραμέτρων $a=1.0,b=4.0,c=1.0,d=1.0$. Η διαφορά σε αυτή την περίπτωση, σε σχέση με την προηγούμενη, είναι το ότι οι δύο σταθερές b και c είναι αρκετά διαφορετικές, δηλαδή η μείωση του πληθυσμού του θηράματος είναι αναλογικά μεγαλύτερη από την αύξηση του πληθυσμού του θηρευτή. Για αυτό το λόγο ο πληθυσμός του θηράματος μεταβάλλεται πολύ περισσότερο από τον πληθυσμό του θηρευτή.

Οι καμπύλες του χώρου φάσεων, ή καλύτερα επιπέδου φάσεων (phase plane, βλέπε επίσης Κεφ. 3), καθώς έχουμε δύο μεταβλητές, δίνονται από την παρακάτω διαφορική εξίσωση:

$$\frac{dv}{du}=\frac{\tilde{d}v(u-1)}{u(1-v)}$$

(6.7)

Ολοκλήρωση στο χώρο φάσεων μας δίνει:

$$\tilde{d}u+v-\ln \mathrm{}(u\tilde{d}v)=H$$

(6.8)

όπου $H$ είναι σταθερά με: $H>H_{\min },H_{\min }=1+\tilde{d}$.

Για $H>H_{\min }$, και συγκεκριμένες τιμές των παραμέτρων του μοντέλου, ο χώρος φάσεων του μοντέλου φαίνεται γραφικά στο Σχήμα 6.2. Κάθε τροχιά (trajectory) είναι περιοδική λύση των αδιάστατων μεταβλητών $u$, και $v$ στο χρόνο $\tau$. Όπως βλέπουμε η μορφή του χώρου φάσεων είναι έλλειψη με κέντρο στο σημείο $(1,1)$.

Τα μοντέλα ανταγωνισμού περιγράφουν ένα οικοσύστημα στο οποίο δύο ή περισσότεροι οργανισμοί ανταγωνίζονται για την εξεύρεση τροφής από την ίδια πηγή.

Έστω δύο είδη με πληθυσμούς $N_1(t)$ και $N_2(t)$ σε χρόνο $t$, τα οποία ανταγωνίζονται για την ίδια πηγή τροφής. Για την δημιουργία του μοντέλου κάνουμε τις παρακάτω υποθέσεις:

• •

  Ο πληθυσμός κάθε είδους, απουσία του άλλου, ακολουθεί λογιστική μορφή.

• •

  Η αλληλεπίδραση (ανταγωνισμός) οδηγεί σε μείωση και των δύο πληθυσμών. Σε πρώτη προσέγγιση, η μείωση αυτή είναι ένας γραμμικός όρος, ως προς το γινόμενο ${\mathrm{{\rm N}}}_1{\mathrm{{\rm N}}}_2$.

Με βάση τις υποθέσεις αυτές, οι εξισώσεις που περιγράφουν την εξέλιξη του πληθυσμού τους είναι:

(6.12)

Στις παραπάνω σχέσεις $r_1,r_2,k_1,k_2,b_{12}$, και $b_{21}$ είναι θετικές σταθερές. Ειδικότερα οι σταθερές $b_{12}$, και $b_{21}$ εκφράζουν το μέτρο του ανταγωνισμού στο συγκεκριμένο σύστημα. Η μελέτη του μοντέλου, όπως και στα μοντέλα θηρευτή-θηράματος, γίνεται σε δύο στάδια.

Πρώτα αδιαστατοποιούμε το σύστημα ορίζοντας τις νέες μεταβλητές:

$$u_1(\tau )=\frac{N_1(t)}{K_1},u_2(\tau )=\frac{N_2(t)}{K_2},\tau =r_{1}t$$

και τις νέες παραμέτροους

$$\rho =\frac{r_2}{r_1},a_{12}=b_{12}\frac{K_2}{K_1},a_{21}=b_{21}\frac{K_1}{K_2}.$$

Με την παραπάνω επιλογή το αρχικό σύστημα εξισώσεων γίνεται:

(6.13)

Τα σημεία ισορροπίας, $u_1^{*}$ και $u_2^{*}$, του παραπάνω συστήματος βρίσκονται λύνοντας τις $f_1=0$ και $f_2=0$ και είναι:

• •

  $(u_1^{*},u_2^{*})=(0,0)$,

• •

  $(u_1^{*},u_2^{*})=(1,0)$,

• •

  $(u_1^{*},u_2^{*})=(0,1)$ και

• •

  $(u_1^{*},u_2^{*})=(\frac{1-a_{12}}{1-a_{12}{a}_{21}},\frac{1-a_{21}}{1-a_{12}{a}_{21}})$.

Η ανάλυση ευστάθειας για κάθε σημείο ισορροπίας γίνεται όπως και στην περίπτωση του μοντέλου Lotka-Volterra, την οποία είδαμε στο προηγούμενο υποκεφάλαιο και αφήνεται ως άσκηση για τον αναγνώστη.

Στο τελευταίο μέρος της μελέτης μοντέλων αλληλεπιδρώντων πληθυσμών εξετάζουμε τα μοντέλα συμβίωσης. Τα μοντέλα αυτά περιγράφουν ένα οικοσύστημα στο οποίο δύο ή περισσότεροι οργανισμοί συμβιώνουν με αρμονικό τρόπο, δηλαδή η αλληλεπίδραση μεταξύ τους αποβαίνει προς όφελος όλων των ειδών. Κλασικό παράδειγμα αποτελεί ένα οικοσύστημα στο οποίο φυτά και πουλιά συμβιώνουν με αρμονικό τρόπο. Έστω και πάλι δύο είδη με πληθυσμούς N1(t) και N2(t) σε χρόνο t, τα οποία συμβιώνουν στο περιβάλλον. Το απλούστερο δυνατό μοντέλο προκύπτει κάτω από τις εξής υποθέσεις:

Υπάρχει άφθονο φαγητό και για τα δύο είδη, δηλαδή Ο πληθυσμός κάθε είδους, απουσία του άλλου, ακολουθεί εκθετική μορφή. Η αλληλεπίδραση (Συμβίωση) οδηγεί σε αύξηση και των δύο πληθυσμών. Σε πρώτη προσέγγιση, η μείωση αυτή είναι ένας γραμμικός όρος, ως προς το γινόμενο Ν1Ν2.

Συνεπώς το μοντέλο το οποίο θα περιέγραφε τους πληθυσμούς θα ήταν:

$$\frac{d{N}_1}{dt}=r_1{N}_1+a_1{N}_1{N}_2,\frac{d{N}_2}{dt}=r_2{N}_2+a_2{N}_1{N}_2$$

(6.14)

όπου $r_1,r_2,a_1$ και $a_2$ είναι θετικές σταθερές.

Το μοντέλο όμως αυτό δεν είναι ρεαλιστικό καθώς προβλέπει συνεχή αύξηση των πληθυσμών $N_1(t)$ και $N_2(t)$, θυμίζει δηλαδή το Μαλθουσιανό-εκθετικό μοντέλο του ενός πληθυσμού.

Ένα πιο ρεαλιστικό μοντέλο συμβίωσης έχει προταθεί [May 1975] σε αναλογία με το μοντέλο ανταγωνισμού που είδαμε στο προηγούμενο υπό-κεφάλαιο. Δηλαδή η διαφορά από το παραπάνω μοντέλο είναι το ότι θεωρούμε πως ο πληθυσμός κάθε είδους, απουσία του άλλου, ακολουθεί λογιστική μορφή. Πιο συγκεκριμένα το μοντέλο αυτό έχει τη μορφή:

(6.15)

όπου $r_1,r_2,k_1,k_2,b_{12}$, και $b_{21}$ είναι θετικές σταθερές.

Ειδικότερα οι σταθερές $b_{12}$, και $b_{21}$ εκφράζουν το μέτρο της συμβίωσης στο συγκεκριμένο σύστημα.

Η αδιαστατοποίηση των παραπάνω εξισώσεων γίνεται με τις νέες μεταβλητές, παρόμοια με το μοντέλο ανταγωνισμού:

$$u_1(\tau )=\frac{N_1(t)}{K_1},u_2(\tau )=\frac{N_2(t)}{K_2},\tau =r_{1}t,\rho =\frac{r_2}{r_1},a_{12}=b_{12}\frac{K_2}{K_1},a_{21}=b_{21}\frac{K_1}{K_2}.$$

Τελικά το σύστημα των αδιάστατων εξισώσεων συμβίωσης γίνεται:

(6.16)

Τα σημεία ισορροπίας, $u_1^{*}$ και $u_2^{*}$, του παραπάνω συστήματος είναι:

• •

  $(u_1^{*},u_2^{*})=(0,0)$,

• •

  $(u_1^{*},u_2^{*})=(1,0)$,

• •

  $(u_1^{*},u_2^{*})=(0,1)$, και

• •

  $(u_1^{*},u_2^{*})=\frac{1+a_{12}}{1-a_{12}{a}_{21}},\frac{1+a_{21}}{1-a_{12}{a}_{21}}$

Η ανάλυση ευστάθειας για κάθε σημείο ισορροπίας αφήνεται ως άσκηση για τον αναγνώστη.