Macierze Lesliego i Markowa  K. Leśniak 1
Model Lesliego
Wyodrębniamy w populacji k grup wiekowych. Po każdej
jednostce czasu następują narodziny i zgony oraz starze-
nie (przechodzenie do następnej grupy wiekowej). Chcemy
symulować zmiany stanu liczebności w poszczególnych gru-
pach wiekowych.
Oznaczmy:
0 mi  liczba potomstwa pojawiającego się co jednost-
kę czasu u osobnika z i-tej grupy wiekowej, i = 1, . . . , k,
[0, 1] si  przeżywalność osobników z i-tej grupy wie-
kowej dożywających przynależności do następnej (i + 1)-ej
grupy wiekowej (jaki procent odobników i-tej grupy dożywa
awansu do (i + 1)-ej grupy wiekowej), i = 1, . . . , k - 1,
xn,i  liczba osobników z i-tej grupy wiekowej w n-tej
jednostce czasu, i = 1, . . . , k.
2 Macierze Lesliego i Markowa  K. Leśniak
Struktura wiekowa populacji w n-tej jednostce czasu
xn = [xn,1, xn,2, . . . , xn,k]T
podlega zależnościom:
k

xn+1,1 = mi � xn,i
i=1
 liczba nowych osobników (wszystkie grupy wiekowe mogą
wydawać potomstwo),
xn+1,i+1 = si � xn,i, i = 1, 2, . . . , (k - 1),
 liczba osobników, które postarzały się awansując do na-
stępnej grupy wiekowej (pozostałe umarły).
Macierzowo: xn+1 = M � xn, gdzie
�ł łł
m1, m2, . . . , mk-1, mk
�ł śł
�ł śł
�ł śł
�ł śł
s1, 0, . . . , 0, 0
�ł śł
�ł śł
�ł śł
�ł śł
M = 0, s2, . . . , 0, 0 .
�ł śł
�ł śł
�ł śł
�ł śł
�ł . . . . . . . . . . . . . . . śł
�ł śł
�ł �ł
0, 0, . . . , sk-1, 0
Powyższy wzór rekurencyjny prowadzi do wzoru  jawnego :
xn = Mn � x0.
Macierze Lesliego i Markowa  K. Leśniak 3
Przykład. W pewnej populacji wyróżniają się dwie grupy
wiekowe:
1  osobniki niedojrzałe, które nie mogą się rozmnażać
(m1 = 0 o przeżywałności s (s1 = s), oraz
2  osobniki dojrzałe o płodności m, które mogą się roz-
mnażać (m2 = m).
�ł łł �ł łł
m1 m2 0 m
�ł śł �ł śł
�ł śł �ł śł
Macierz Lesliego: M = = .
�ł �ł �ł �ł
s1 0 s 0
Przez indukcję można zauważyć, że:
�ł łł �ł łł
0 mn+1 � sn mn � sn 0
�ł śł �ł śł
�ł śł �ł śł
M2n+1 = , M2n =
�ł �ł �ł �ł
sn+1 � mn 0 0 sn � mn .
Zatem
�ł łł �ł łł
0 mn+1 � sn x0,1
�ł śł �ł śł
�ł śł �ł śł
x2n+1 = M2n+1 � x0 = �
�ł �ł �ł �ł
sn+1 � mn 0 x0,2 =
�ł łł
mn+1 � sn � x0,2
�ł śł
�ł śł
= �
�ł �ł
sn+1 � mn � x0,1 = (m � s)n [m � x0,2, s � x0,1]T,
�ł łł �ł łł
mn � sn 0 x0,1
�ł śł �ł śł
�ł śł �ł śł
x2n = M2n �x0 = �x0.
�ł �ł �ł �ł
0 sn � mn � x0,2 = (m�s)n
Zachowanie asymptotyczne:
�ł łł ńł
�ł
0
�ł
�ł
xn,1 [ ] , gdy 0 < ms < 1,
�ł śł
0
�ł śł
xn =
�ł �ł -
�ł
"
�ł
ół
xn,2 n" [ ] , gdy ms > 1.
"
4 Macierze Lesliego i Markowa  K. Leśniak
Przykład (cykliczne zmiany struktury wieku).
�ł łł
0 0 4
�ł śł

�ł śł
1 0 0
�ł śł
1
�ł śł
M = 0 0 . A
atwo dostrzec, że: M3 = 0 1 0 .
�ł śł
�ł 3 śł
0 0 1
�ł �ł
3
0 0
4
Zachowanie asymptotyczne:
x3n = x0 = [x01, x02, x03]T,
�ł łł
T
1 3
�ł śł
�ł �ł
x3n+1 = x1 = [x11, x12, x13]T = 4 x03, x01, x02 ,
3 4
�ł łł
T
4 1
�ł śł
�ł �ł
x3n+2 = x2 = [x21, x22, x23]T = 3 x02, x03, x01 .
3 4
Cykl powtarza się co trzy jednostki czasu:
x0, x1, x2, x0, x1, . . .
Macierze Lesliego i Markowa  K. Leśniak 5
Przykład (stabilna struktura wieku).
�ł łł
2 1 1
�ł śł
�ł śł
�ł śł
�ł śł
M = 0.25 0 0 .
�ł śł
�ł śł
�ł �ł
0 0.(3) 0
W tabeli zbieramy dane o procentowym udziale w liczeb-
ności populacji poszczególnych grup wiekowych po trzech
iteracjach.
udział procentowy
liczebność

3
[xn,i/ xn,i, i = 1, 2, 3]
początkowa
i=1
[1, 0, 0] [49, 772%, 26, 941%, 23, 288%]
[0, 20, 0] [49, 514%, 27, 185%, 23, 301%]
[0, 0, 8] [50, 000%, 24.999%, 24, 999%]
[100, 100, 100] [49, 749%, 26, 936%, 23, 318%]
[250, 10000, 100] [49, 561%, 27, 138%, 23, 300%]
[100000, 5000, 10] [49, 771%, 26, 942%, 23, 288%]
[200, 45000, 6000] [49, 534%, 27, 133%, 23, 334%]
Struktura wiekowa populacji w stosunkowo krótkim czasie
stabilizuje się i zależy od takich parametrów specyficznych
dla populacji jak śmiertelność i rozrodczość.
6 Macierze Lesliego i Markowa  K. Leśniak
Macierze Markowa
�ł łł
p11 p12 p13 . . . p1s
�ł śł
�ł śł
�ł śł
�ł śł
p21 p22 p23 . . . p2s
�ł śł
�ł śł
M = = [pij]1 i,j s
�ł śł
�ł śł
�ł . . . . . . . . . . . . . . . śł
�ł śł
�ł �ł
ps1 ps2 ps3 . . . pss
 macierz stochastyczna (= macierz Markowa
= macierz przejścia), gdy
s

"i,j 0 pij 1, "i pij = 1.
j=1
Terminologia:
" {1, . . . , s}  zbiór stanów układu;
" pij  prawdopodobieństwo przejścia
ze stanu i w stan j, i, j " {1, . . . , s};
" M = [pij]i,j"{1,...,s}  macierz przejścia
dla jednorodnego łańcucha Markowa;
" p(n)  prawdopodobieśtwo z jakim cząstka, która jest w
ij
stanie i znajdzie się w stanie j po n cyklach.
Macierze Lesliego i Markowa  K. Leśniak 7
Przykład. Niech (1) czerwony i (2) niebieski będą stana-
mi jakie przyjmuje pewna cząstka. Zmienia ona swój stan
zgodnie z macierzą przejścia
�ł łł
1 2
�ł śł
�ł śł
3 3
�ł śł
M = .
�ł śł
�ł �ł
3 1
4 4
W szczególności:
2
p12 =  prawdopodobieństwo zmiany
3
z czerwonego na niebieski,
3
p21 =  prawdopodobieństwo zmiany
4
z niebieskiego na czerwony.
(a) Wyznaczyć prawdopodobieństwo p(2) znalezienia się
11
czerwonej cząstki w stanie czerwoności po dwóch cyklach.

1


p11 p12







p(2) 2 I cykl
11
1




p11 p21






II cykl
1

1 1 2 3 11 1
p(2) = p11 � p11 + p12 � p21 = � + � = = p11.
11
3 3 3 4 18 3
8 Macierze Lesliego i Markowa  K. Leśniak
(b) Wyznaczyć wszystkie prawdopodobieństwa p(3) prze-
ij
miany ze stanu i w stan j po trzech cyklach.

i


pi1 pi2







p(2) 2 I cykl
ij
1




p1j p2j






II cykl
j

2

p(2) = pi1 � p1j + pi2 � p2j = pik � pkj,
ij
k=1
macierz przejścia po dwóch cyklach:
[p(2)]i,j = M � M = M2.
ij
Macierze Lesliego i Markowa  K. Leśniak 9

i


p(2) p(2)

i1 i2






p(3) 2 II cykl
ij
1




p1j p2j






III cykl
j

2

p(3) = p(2) � p1j + p(2) � p2j = p(2) � pkj,
ij i1 i2
ik
k=1
macierz przejścia po trzech cyklach:
[p(3)]i,j = [p(2)]i,k � [pkj]k,j = M2 � M = M3.
ij
ik
10 Macierze Lesliego i Markowa  K. Leśniak
(c) Prawdopodobieństwo pozostawania czerwonym przez
n cykli wynosi
n
p11 � p11 � . . . � p11 = p11 .

n razy
Prawdopodobieństwo pozostawania niezmiennie czerwonym
wynosi
�ł �ł
n
1
n �ł �ł
�ł łł
lim p11 = lim = 0.
n" n"
3
Co dzieje się z prawdopodobieństwem p(n) powrotu do
11
stanu czerwoności po n cyklach?

i


p(n-1) p(n-1)

i1 i2






p(n) 2 cykl
ij
1

(n - 1)



p1j p2j






cykl n
j

2

p(n) = p(n-1) � p1j + p(n-1) � p2j = p(n-1) � pkj,
ij i1 i2
ik
k=1
macierz przejścia po n cyklach:
[p(n)]i,j = [p(n-1)]i,k � [pkj]k,j = Mn-1 � M = Mn.
ij
ik
Macierze Lesliego i Markowa  K. Leśniak 11
Aby znalez
ć p(n) musimy wyznaczyć potęgę Mn.
11
�ł łł
ąn �n
�ł śł
�ł śł
Oznaczmy: Mn =
�ł �ł
łn �n . Wówczas
�ł łł �ł łł �ł łł
ąn+1 �n+1 ąn �n ą1 �1
�ł śł �ł śł �ł śł
�ł śł �ł śł �ł śł
�ł �ł �ł �ł �ł �ł
łn+1 �n+1 = łn �n � ł1 �1 ,
skąd
ńł �ł
�ł �ł
�ł �ł
�ł ąn + �n = 1 (1)
�ł
�ł �ł
�ł �ł
�ł żł
�ł
�ł
�ł
�ł
(*)
�ł
�ł �ł
1 3
�ł �ł
�ł �ł
�ł ąn+1 = ąn + �n (2) �ł
�ł �ł
�ł �ł
3 4
�ł
�ł
�ł
�ł
�ł
�ł
2 1
�ł
�ł
�ł
�n+1 = ąn + �n
�ł
�ł
3 4
�ł
�ł
�ł
�ł
�ł
�ł
�ł
1 2
�ł
�ł
�ł ą1 = , �1 = (3)
�ł
3 3
�ł
�ł
�ł
�ł
�ł
łn + �n = 1
�ł
�ł
�ł
�ł
�ł
�ł
�ł
�ł
�ł
1 3
�ł
�ł
�ł łn+1 = łn + �n
�ł
�ł
3 4
�ł
�ł
�ł
�ł
�ł
�ł
2 1
�ł
�ł
�ł
�n+1 = łn + �n
�ł
�ł
3 4
�ł
�ł
�ł
�ł
�ł
�ł
�ł
3 1
�ł
ół
ł1 = , �1 = .
4 4
Wystarczy rozwiązać podukład (*). Podstawiając (1) do
(2) dostajemy równanie(L-1):
1 3 5 3
ąn+1 = ąn + (1 - ąn) = - ąn + ,
3 4 12 4
którego rozwiązanie ogólnym jest
�ł �ł
n
5
�ł �ł
�ł �ł
�ł łł
1 - -
n
5 3
�ł �ł 12
�ł łł �ł �ł
ąn = - � ą0 + � .
5
12 4
�ł �ł
�ł łł
1 - -
12
12 Macierze Lesliego i Markowa  K. Leśniak
Uwzględniając warunek początkowy (3):
�ł �ł
1 5 3
�ł �ł
�ł łł
= ą1 = - ą0 + �! ą0 = 1
3 12 4
dostajemy
�ł �ł
n
5
�ł �ł
�ł �ł
�ł łł
1 - -
n
5 3
�ł �ł 12
�ł łł �ł �ł
ąn = - + � .
5
12 4
�ł �ł
�ł łł
1 - -
12
Zatem
3
9
p(n) = ąn n" 4 = .
-
11
5
17
1 +
12
Podobnie
8
p(n) = �n = 1 - ąn n" .
-
12
17
Przyjmując za (ąn, �n) odpowiednio ciągi (łn, �n) otrzymu-
jemy rozwiązanie ogólne tego samego kształtu: �n = 1 - łn,
�ł �ł
n
5
�ł �ł
�ł �ł
�ł łł
1 - -
n
5 3
�ł �ł 12
�ł łł �ł �ł
łn = - � ł0 + � .
5
12 4
�ł �ł
�ł łł
1 - -
12
Po uwzględnieniu warunków początkowych
�ł �ł
n
5
�ł �ł
�ł �ł
�ł łł
1 - -
n
3 5 3
�ł �ł 12
�ł łł �ł �ł
= ł1 = - � ł0 + � �! ł0 = 0
5
4 12 4
�ł �ł
�ł łł
1 - -
12
Macierze Lesliego i Markowa  K. Leśniak 13
uzyskujemy
�ł �ł
n
5
�ł �ł
�ł łł
1 - -
3
12
�ł �ł
łn = � .
5
4
�ł �ł
�ł łł
1 - -
12
Zatem jak poprzednio
3
9
p(n) = łn n" 4 = ,
-
21
5
17
1 +
12
8
p(n) = �n n" .
-
22
17
W ten sposób otrzymaliśmy opis asymptotycznego zacho-
wania się macierzy przejścia po n cyklach:
�ł łł
9 8
�ł śł
�ł śł
17 17
�ł śł
Mn - .
�ł śł
�ł �ł
n"
9 8
17 17
Mówimy, że macierz M jest ergodyczna, a prawdopodo-
bieństwa graniczne
9
lim p(n) = lim p(n) = ,
11 21
17
n" n"
8
lim p(n) = lim p(n) = ,
12 22
17
n" n"
nazywane są prawdopodobieństwami ergodycznymi.
14 Macierze Lesliego i Markowa  K. Leśniak
(d) Niech Ą = [Ą1, Ą2] będzie początkowym rozkładem
prawdopodobienstwa: zanim zadziałał mechanizm zmian ko-
lorów cząstka znajdowała się w stanie (1) z prawdopodo-
bieństwem Ą1, a w stanie (2) z prawdopodobieństwem Ą2.
Wyznaczyć rozkład prawdopodobieństwa stanu cząstki Ą(1) =
(1) (1)
[Ą1 , Ą2 ] po jednym cyklu.


Ą1 Ą2







stan
1 2

początkowy

p11 p12 p21 p22





stan

po 1 cyklu
1 2 1 2

(1)
Ą1 = Ą1 � p11 + Ą2 � p21,
(1)
Ą2 = Ą1 � p12 + Ą2 � p22,
czyli
�ł łł
p11 p12
(1) (1) �ł śł
�ł śł
[Ą1 , Ą2 ] = [Ą1, Ą2] �
�ł �ł
p21 p22 = [Ą1, Ą2] � M.
Macierze Lesliego i Markowa  K. Leśniak 15
(e) Przez indukcję rozkład prawdopodobieństwa Ą(n) =
(n) (n)
[Ą1 , Ą2 ] stanu cząstki po n cyklach opisuje zależność
�ł łł
p11 p12 n
(n) (n) �ł śł
�ł śł
[Ą1 , Ą2 ] = [Ą1, Ą2] �
�ł �ł
p21 p22 = [Ą1, Ą2] � Mn.
Uzasadnienie: ćwiczenie.
Niezależnie od rozkładu początkowego Ą
Ą(n) = Ą � Mn - Ą",
n"
" "
gdzie Ą" = [Ą1, Ą2],
ńł
�ł
9
�ł
�ł
�ł , j = 1
�ł
17
"
Ąj = lim p(n) =
ij �ł
�ł
n"
8
�ł
�ł
ół
, j = 2
17
niezależne od i = 1, 2 prawdopodobieństwa graniczne zna-
lezione w punkcie (c). Tym samym rozkład graniczny Ą" wy-
stępuje we wszystkich wierszach macierzy granicznej lim Mn.
n"
16 Macierze Lesliego i Markowa  K. Leśniak
Przykład.
�ł łł
0.2 0.4 0.4
�ł śł
�ł śł
�ł śł
�ł śł
M = 0.25 0.25 0.5 .
�ł śł
�ł śł
�ł �ł
0.(3) 0.(3) 0.(3)
W tabeli zbieramy dane z czterech kroków iteracji wyko-
nanych dla trzech rozkładów początkowych.
n Ą(n) Ą(n) Ą(n)
0 [0, 0, 1] [0.5, 0.25, 0.25] [.667, .333, 0]
1 [.333, .333, .333] [.246, .346, .408] [.217, .350, .433]
2 [.261, .327, .411] [.272, .321, .407] [.275, .318, .406]
3 [.271, .323, .405] [.271, .325, .406] [.270, .325, .404]
4 [.270, .324, .405] [.271, .324, .406] [.270, .324, .406]
Ą" = [0.(270), 0.(324), 0.(405)].
Macierze Lesliego i Markowa  K. Leśniak 17
Sam rozkład ergodyczny Ą" łatwo znalez
ć nie iterując ani
macierzy ani rozkładów, gdyż jest on stacjonarny.
Twierdzenie 1
(stacjonarność rozkładu ergodycznego)
Jeżeli łańcuch Markowa M jest ergodyczny tzn.
n-te rozkłady prawdopodobieństwa Ąn = Ą Mn zbiegają
do pewnego rozkładu Ą" niezależnie od rozkładu począt-
kowego Ą, to
Ą" = Ą" � M
jest jedynym rozkładem stacjonarnym.
Dowód. Stacjonarność:
Ą(n) = Ą � Mn - Ą",
n"
Ą" �!- Ą(n+1) = Ą Mn+1 = (Ą Mn) M - Ą" M.
n" n"
Jedyność:
Ą = Ą M �!
Ć Ć
Ą = Ą M = (Ą M) M = Ą M2 = . . . = Ą Mn - Ą"
Ć Ć Ć Ć Ć
n"
Ć
�! Ą = Ą".

18 Macierze Lesliego i Markowa  K. Leśniak
�ł łł
1 2
�ł śł
�ł śł
3 3
�ł śł
Przykład. M = , Ą = [Ą1, Ą2], Ąi 0, Ąi = 1.
�ł śł
�ł �ł
3 1
i
4 4
�ł łł
1 3 2 1
�ł śł
�ł �ł
[Ą1, Ą2] = [Ą1, Ą2] M = Ą1 + Ą2, Ą1 + Ą2
3 4 3 4
9 8
�! Ą1 = , Ą2 = 1 - Ą1 = .
17 17
Jednak nawet jedyność rozkładu stacjonarnego nie gwa-
rantuje jego ergodyczności.
Przykład (łańcuch okresowy).
�ł łł
0 1
�ł śł

�ł śł
M = , Ą = [Ą1, Ą2], Ąi 0, Ąi = 1.
�ł �ł
i
1 0
[Ą1, Ą2] = [Ą1, Ą2] M = [Ą2, Ą1]
1
�! Ą1 = Ą2 = .
2
Znaleziony rozkład nie jest ergodyczny, bo
ńł �ł łł
�ł
�ł
�ł
1 0
�ł �ł śł
�ł
�ł �ł śł
�ł
, n  parzyste,
�ł �ł
Mn = 0 1
�ł
�ł
�ł
�ł
�ł
�ł
ół
M, n  nieparzyste.
Jak się przekonać, że znaleziony rozkład stacjonarny jest
ergodyczny?
Macierze Lesliego i Markowa  K. Leśniak 19
Jako że znajdowanie jawnej postaci n-tej potęgi macie-
rzy jest pracochłonne (nawet z użyciem postaci Jordana),
powstaje pytanie, czy nie można określić ergodyczności ma-
cierzy Markowa bezpośrednio na podstawie jej współczyn-
ników.
Takiego kryterium dostarcza
Twierdzenie 2 (ergodyczne dla łańcuchów)
Niech macierz Markowa M = [pij]i,j"{1,...,s}
spełnia "i,j pij > 0.
Wówczas
" "
"j " Ąj > 0 "i n" p(n) = Ąj ,
lim
ij

s

" "
Ąj = 1, Ą" = Ąj s = Ą" M.
j=1
j=1

Przypomnijmy, że Mn = p(n) , a więc
ij
i,j"{1,...,s}
Mn - [Ą", . . . , Ą"]T,
n"
s razy
Ą(n) = Ą Mn - Ą".
n"
20 Macierze Lesliego i Markowa  K. Leśniak
Dowód. Ustalmy 0 < � < min pij
i,j
i przyjmijmy oznaczenia:
[!]
(n)
ąj = min p(n) > 0,
ij
i=1,...,s
(n) (n)
�j = max p(n) ąj .
ij
i=1,...,s
[!] Jeśli "ij pij > 0, to "n "ij p(n) > 0 (ćwiczenie).
ij
Oczywiście
[p(n+1)]i,j = Mn+1 = Mn � M = [pil]i,l � [p(n)]l,j,
ij
lj

p(n+1) = pil � p(n).
ij
lj
l
Stąd


p(n+1) = pil - � � p(n) � p(n) + � � p(n) � p(n)
ij
jl lj jl lj
l l


(n)
�j � pil - � � p(n) + � � p(2n) =
jj
jl
l
�ł �ł
�ł �ł
�ł �ł
�ł �ł

(n)
�ł �ł
�ł �ł
= �j � pil - � � p(n) + � � p(2n),
�ł �ł jj
jl
�ł �ł
l l
�ł łł

=1 =1
czyli
(n+1) (n)
�j = max p(n+1) �j � (1 - �) + � � p(2n).
ij jj
i
Analogicznie
(n+1) (n)
ąj ąj � (1 - �) + � � p(2n) (*).
jj
A
ącznie:

(n+1) (n+1) (n) (n)
�j - ąj (1 - �) � �j - ąj ,
Macierze Lesliego i Markowa  K. Leśniak 21
skąd

(n) (n) (1) (1)
�j - ąj (1 - �)n-1 � �j - ąj n" 0.
-
Przechodząc w (*) z � 0 widzimy, że
(n) (n+1)
ąj ąj 1,

"
(n)
tzn. ąj n=1  rosnący i ograniczony. Istnieją więc granice
(n)
"
Ąj = lim ąj ,
n"
"
przy czym Ąj > 0, bo
(n) (n) (1)
lim ąj ąj ąj = min pij min pij > � > 0.
n"
i ij
Ostatecznie


(n) (n)
"


p(n) - Ąj �j - ąj n" 0,
-
ij
"
a tym samym istnieją lim p(n) = Ąj > 0 oraz
ij
n"

"
Ąj = lim p(n) = lim p(n) = 1.
ij ij
n" n"
j j j

Uwaga: Wystarczy założyć, że dla pewnego n0 potęga
Mn0 ma wszystkie wyrazy dodatnie: p(n0) > 0.
ij