Evrimsel oyun teorisi ve seçilim

Martin Nowak evolutionary dynamics kitabında üremeyi şöyle basitçe modeller:

\(\dot{x} = \frac{dx}{dt} = r\cdot x \)

Burada \(x \) değişkeni türün nüfusunu belirtiyor. Dolayısıyla \(\dot{x} \) terimi de türün nüfusunun zamana göre değişimini gösteriyor, yani zaman geçtikçe nüfusunun ne kadar arttığını ya da azaldığını ifade ediyor. Bu artış ya da azalış da o anki nüfusun sabit bir \(r \) sayısıyla çarpımına eşit olduğu bir modelle gösterilsin diyor. Bu arada doğum ve ölüm oranlarının çoğalma katsayısında olduğunu da unutmamak lazım. Yani aradaki fark çoğalmayı veriyor ve bu fark aslında nüfusun azalması da demek olabilir.

Bu şekildeki bir modelde nüfus artışı tabii ki üssel bir davranış gösterecektir.

\(r = 1.1 \) değeri için yukarıdaki grafiği elde ettim. Toplamda yirmi adım bile gidilmeden nüfus bir milyonu geçti. Böyle bir ortamda bütün türler aynı hızda çoğalırdı ama. Nowak hemen sonra gelen seçilim kısmının başında “seçilim değişik türlerin bireyli değişik oranlarda çoğalırsa meydana gelir” diyor. Eğer her tür aynı oranda çoğalsaydı nüfuslarını birbirlerinden ayıran tek fark başlangıçtaki nüfusları olacaktı. Diğer bir deyişle nüfusların başlangıçtaki oranları nesiller geçtikçe hiç değişmeyecekti. Bu durumda gerçekten de bir seçilimden bahsetmek zor olabilir.

İkinci cümlede ise seçilim olması için en az iki tür olmalıdır diyor. Bu da seçimden bahsedilebilmesi için anlaşılır bir varsayım hatta tanım bile olabilir.

O zaman iki türden oluşan modeli şöyle veriyor:

\(\dot{x} = a\cdot {x}\)

\(\dot{y} = b\cdot {y}\)

Burada \(a > b \) ise birinci tür ikinci türden çok daha hızlı üreyecektir ve zaman geçtikçe birinci türün nüfusunun ikinci türün nüfusuna oranı sürekli artacaktır. \(b > a \) olduğu durumda da tersi senaryo gözlenecektir.

Tabii ki her türün sınırsız büyüdüğü modeller pek gerçekçi değil. Ekosistemimizin maksimum birey sayısına sahip olduğunu varsayalım ve modelimizi buna uygun hale getirelim. Bu modelde toplam bir sayıdan bahsetmek yerine nüfus oranları kullanılacak. Yine aynı değişkenleri kullanıyoruz ama artık \(x \) değişkeni birinci türün nüfusunu değil de birinci türün nüfusunun toplam birey sayısına oranını veriyor. Bu oranlar haliyle minimum 0 ve maksimum 1 değerine sahip olabilirler.

Kitapta bunun için şu model veriliyor:

\(\dot{x} = a\cdot {(x-\phi)}\)

\(\dot{y} = b\cdot {(y- \phi)}\)

Buradaki \(\phi \) teriminin görevi türlerin nüfuslarının oranlarının toplamını 1 değerinde sabit tutmak. Bunu sağlamak için de

\(\phi = a\cdot {x} + b\cdot{y} \)

eşitliğinin sağlanması gerekiyor. Yani her adımda \(\phi \) değeri o adımdaki nüfusların oranına göre tekrar hesaplanıyor.

Aşağıdaki python programıyla bu modeli denedim.

import numpy as np
import matplotlib.pyplot as plt


a = 1.1
b = 1.3
number_of_iterations = 100

def constrainedA(x, phi) :
    return x*(a - phi)

def constrainedB(y, phi) :
    return y*(b - phi)

def phi(x, y):
    return a*x + b*y;
x = 0.5
y = 0.5

x_population = np.array([[0, x]])
y_population = np.array([[0, y]])
total = np.array([[0, x+y]])

for i in range(0, number_of_iterations):
    phi_value = phi(x, y)
    delta_x = constrainedA(x, phi_value)
    delta_y = constrainedB(y, phi_value)
    x += delta_x
    y += delta_y
    x_population = np.append(x_population, np.array([[i, x]]), axis = 0)
    y_population = np.append(y_population, np.array([[i, y]]), axis = 0)
    total = np.append(total, np.array([[i, x+y]]), axis = 0)

plt.plot(x_population[:, 0], x_population[:, 1])
plt.plot(y_population[:, 0], y_population[:, 1])
plt.plot(total[:, 0], total[:, 1])
plt.xlabel("zaman")
plt.ylabel("nüfus")
#plt.plot(y_population[:, 0], y_population[:, 1])
plt.show()

Sonuç olarak da şu grafiği elde ettim:

Bu programda adım sayısını 20’den 100’e yükselttim. Bu sayede türlerin birinin diğerini ekosistemden nasıl sildiğini görebiliyoruz. Kırmızı çizgiyle gösterilen türün çoğalma katsayısı 1.3, mavi çizgiyle gösterilen türünkü ise 1.1 idi. Yeşi çizgi de nüfus oranlarının toplamını gösteriyor ve sürekli 1 değerine sahip.

Martin Nowak bu modeli daha iyinin hayatta kalmasına (Survival of the fitter) örnek olarak veriyor. Gerçekten de bu modelde daha iyi çoğalma yeteneğine kalan tür hayatta kalıyor ve diğeri yok oluyor. Sonra aynı modeli ikiden fazla tür için kurup bu sefer de en iyinin hayatta kalması (Survival of the fittest) fikrini gösteriyor. Bu modelde iki genellemeye gidiyor:

\(\phi = \sum_{i=1}^{n} x_{i} \cdot {f_{i}} \)

\(\dot{x_{i}} = x_{i} \cdot (f_{i} – \phi) \)

Bu modelde \(x_{i}\) i numaralı türün nüfus oranını \(f_{i}\) de aynı türün çoğalma fonksiyonunu (fitness) gösteriyor. Modelin gerisi aynı şekilde çalıştırılıyor.

Bu noktadan sonra Nowak iki yeni durum için modelimiz nasıl olmalı diye bir soru soruyor:

  1. Ekosistemdeki ilk mevcut tür daha sonradan gelen tür ne kadar iyi olursa olsun mücadeleyi kazansın (Survival of the first)
  2. Ekosistemdeki her tür hayatta kalsın (Survival of all)

Bu durumda Nowak doğrusal fitness fonksiyonu şartından vazgeçmemiz gerekir diyor ve şu modeli sunuyor:

\(\dot{x} = a\cdot {x^c}-\phi \cdot {x}\)

\(\dot{y} = b\cdot {y^c}- \phi \cdot {y}\)

Bu modelde nüfus oranlarının toplamını sabitlemek için şu eşitliğe ihtiyacımız var:

\(\phi = a\cdot{x^c} + b\cdot {y^c} \)

Modeldeki \(c \) sabitinin birden büyük seçersek 1. olasılığı modellemiş oluyoruz. 1. senaryonun aklıma takılan kısmı bunu nasıl programlayacağımdı. Bu senaryoda söylenen şey şu: Mesela bütün toplam nüfus sadece \(x \) türünden oluşuyorsa sisteme ekleyeceğimiz bir \(y \) türü bireyi çoğalma potansiyeli diğer türden yüksek olsa bile yok olmaya mahkumdur. Yani programda ilk anda \(x = 1 \) ile başlamam gerekecek. Ardından sisteme minimum bir \(y \) değeri eklemem lazım. Bunun yerine başlangıçta minimum \(y \) ile başlamayı seçtim.

import numpy as np
import matplotlib.pyplot as plt


a = 1.1
b = 2.0
c = 1.3
number_of_iterations = 20

def constrainedA(x, phi) :
    return a*pow(x, c) - x*phi

def constrainedB(y, phi) :
    return b*pow(y, c) - y*phi

def phi(x, y):
    return a*pow(x, c) + b*pow(y, c);
x = 0.999
y = 0.001

x_population = np.array([[0, x]])
y_population = np.array([[0, y]])

for i in range(0, number_of_iterations):
    phi_value = phi(x, y)
    delta_x = constrainedA(x, phi_value)
    delta_y = constrainedB(y, phi_value)
    x += delta_x
    y += delta_y
    x_population = np.append(x_population, np.array([[i, x]]), axis = 0)
    y_population = np.append(y_population, np.array([[i, y]]), axis = 0)

plt.plot(x_population[:, 0], x_population[:, 1])
plt.plot(y_population[:, 0], y_population[:, 1])
plt.show()

Bu sefer programda toplam nüfus oranlarını göstermek istemedim, çünkü o zaman birinci türün çok çabuk 1 değerine ulaşması net gözükmüyordu.

Grafikte çok net görülmüyor ama daha ilk adımda birinci tür bütün popülasyonu ele geçirdi.

Bu senaryo tabii ki her başlangıç değeri için bu sonucu vermiyor. Örneğin aynı \(c \) sabiti için \(x \) türünü nüfusun yüzde sekseni olacak şekilde başlatınca bu sonuç çıkıyor.

Görüldüğü gibi daha iyi fitness fonksiyonuna sahip olan \(y \) yok olmadığı gibi bütün ekosistemi de ele geçirdi.

Eğer \(c < 1 \) olacak şekilde bir seçim yaparsak iki tür de yok olmadan beraber yaşama şansını yakalayabiliyor.

Bunun denemesini de aşağıdaki programla yaptım.

import numpy as np
import matplotlib.pyplot as plt


a = 1.2
b = 0.8
c = 0.5
number_of_iterations = 20

def constrainedA(x, phi) :
    return a*pow(x, c) - x*phi

def constrainedB(y, phi) :
    return b*pow(y, c) - y*phi

def phi(x, y):
    return a*pow(x, c) + b*pow(y, c);
x = 0.5
y = 0.5

x_population = np.array([[0, x]])
y_population = np.array([[0, y]])

for i in range(0, number_of_iterations):
    phi_value = phi(x, y)
    delta_x = constrainedA(x, phi_value)
    delta_y = constrainedB(y, phi_value)
    x += delta_x
    y += delta_y
    x_population = np.append(x_population, np.array([[i, x]]), axis = 0)
    y_population = np.append(y_population, np.array([[i, y]]), axis = 0)

plt.plot(x_population[:, 0], x_population[:, 1], label='x')
plt.plot(y_population[:, 0], y_population[:, 1], label='y')
plt.legend()
plt.show()

Grafikte de görüldüğü gibi iki tür de yok olmadı.

Heralde bu özellikleri sağlayan daha başka modeller de vardır ama konuya bu şekilde girilmesi hoşuma gitti. Bir modelin olası sonuçlarından çok, istenen sonucu açıklayabilecek bir model sunma tekniği ilgimi çekti. Çalışmalarıma biraz da bu kitaptan devam edeyim, belki şansım bu sefer döner.

Network karşılıklılığı (5)

Bu tür oyunlarda işbirlikçiler duruma göre ihanetçilerle beraber yaşayabiliyor ve hatta bazı durumlarda daha üstün bile olabiliyorlar. Anladığım kadarıyla burada Oyuncular bir çizgenin (graph) köşeleri olacak şekilde tanımlanıyor. Bu durumda her köşe ya işbirlikçi ya da ihanetçi stratejilerine sahip oluyor. Sonra bir \(k \) tamsayısı ile her köşenin kaç kenara sahip olduğu tanımlanıyor. Yani her oyuncu bu kadar komşuya sahip oluyor. Her oyuncu sadece kendi komşularıyla oynuyor ve standard kazanç matrisine göre kazançlar hesaplanıyor. Bir oyuncu birden fazla komşuyla oynadığında kazançları da oyunlardaki kazançların toplamı oluyor.

Bu modelin evrimleşmesi de şu şekilde gerçekleşiyor. Herhangi bir adımda seçilen rastgele bir oyuncu ölüyor. Sonra bu ölen oyuncunun köşesini kapmak için o köşenin komşu oyuncularından birisi kazançlarına göre seçiliyor.

Her komşu iki oyuncu birbiriyle şu matrise göre oyun oynuyor:

\(M = \begin{bmatrix} R&&S \\T&&P \end{bmatrix} \)

Yukarıdaki kurallarla birleştirince matris şu biçime dönüşüyor:

\(M = \begin{bmatrix} R&&S + H\\T – H&&P \end{bmatrix} \)

\(H = \frac{(k+1)(R – P) – T + S}{(k+1)(k – 2)} \)

Aşağıdaki python programıyla bu oyunun simülasyonunu yapmaya çalıştım:

import numpy as np
import matplotlib.pyplot as plt

strategies = np.array([[1, 0], [0, 1]])

k = 3 #number of edges

R = 3
S = 5
T = 5
P = 1

H = ((k + 1) * (R - P) - T + S)/((k + 1)*(k - 2))
payoff = np.array([[R, S + H], [T - H, P]])

number_of_iterations = 100

increment = 0.01
steady_states = np.array([[0, 0]])


for s1 in np.arange(0, 1.0, increment):
    s2 = 1 - s1
    species = np.array([s1, s2])

    for i in range(0, number_of_iterations):

        difference = (payoff.dot(species) - species.dot(payoff).dot(species))*species
        species = species + difference
        species = np.clip(species, 0, 1)

        
    steady_states = np.append(steady_states, np.array([[s1, species[0]]]), axis = 0)

plt.plot(steady_states[:, 0], steady_states[:, 1])
plt.ylabel("işbirlikçi türün sondaki oranı")
plt.xlabel("işbirlikçi türün başlangıçtaki oranı")
plt.show()

Bu simülasyonun sonucunda da aşağıdaki grafiği elde ettim.

k = 3 için elde ettiğim sonuç

Bu kazanç matrisi ve model için işbirlikçilik daha iyi bir stratejiymiş.

Dolaylı karşılıklılık

Doğrudan karşılıklılıkta birisine karlı vereceğimiz kararı, o kişinin bize geçmişte nasıl davrandığına bakarak veriyorduk. Dolaylı karşılıklılıkta bu kararı o birisinin geçmişte başkalarına nasıl davrandığına bakarak veriyoruz.

Oyuncular arasında tutsak ikilemindeki gibi oyunlar oynanıyor. Bu karşılaşmalar da diğer oyuncular tarafından izleniyor. Karşılaşmalarda işbirliğinin getirisi az. Buna karşın ihanetin kısa vadede getirisi fazlayken aynı zamanda kötü ün de kazandırıyor.

Karşılaşmaların oynanması şu şekilde. İhanet edenler kimle oynarsa oynasın ihanet ediyor. İşbirlikçiler işbirlikçilere karşı her zaman işbirliği yapıyor. İşbirlikçiler bir ihanetçiyi \(q \) ihtimalle tanıyor ve ihanetçiyle sadece \(1 – q \) ihtimalle işbirliği yapıyor. Bu şartlar altında payoff matrisi aşağıdaki şekle dönüşüyor.

\(M = \begin{bmatrix} R&&(1-q)S+qP\\(1-q)T+qP&&P \end{bmatrix} \)

Bir başkasının ününü bilme şansı \(q \) aşağıdaki eşitsizliği sağladığında doğrudan karşılıklılıkla aynı sonuçlar elde ediliyor.

\(q > \frac{T-R}{T-P} \)

Makaledeki bu sonuçları denemek için aşağıdaki python programını yazdım ve denemeler yaptım.

import numpy as np
import matplotlib.pyplot as plt

strategies = np.array([[1, 0], [0, 1]])

R = 5
S = 1
T = 7
P = 3
q = 0.51
prisoner_dilemma_payoff = np.array([[R, (1 - q)*S + q*P], [(1-q)*T + q*P, P]])

number_of_iterations = 100

payoff = prisoner_dilemma_payoff

increment = 0.01
steady_states = np.array([[0, 0]])


for s1 in np.arange(0, 1.0, increment):
    s2 = 1 - s1
    species = np.array([s1, s2])

    for i in range(0, number_of_iterations):

        difference = (strategies.dot(payoff).dot(species) - species.dot(payoff).dot(species))*species
        species = species + difference
        species = np.clip(species, 0, 1)
        
    steady_states = np.append(steady_states, np.array([[s1, species[0]]]), axis = 0)

plt.plot(steady_states[:, 0], steady_states[:, 1])
plt.ylabel("1. türün sondaki oranı")
plt.xlabel("1. türün başlangıçtaki oranı")
plt.show()

Programdaki R, S, T ve P değerleri için yukarıdaki eşitsizlik aşağıdaki gibi oluyor.

\(q > \frac{7 – 5}{7 – 3} = \frac {2}{4} = 0.5 \)

Yani programdaki \(q = 0.51\) değeri için işbirlikçi stratejinin kararlı olduğu bir dağılım da olmasını bekliyordum. Programı çalıştırınca aşağıdaki sonucu aldım.

\(q = 0.51\)

Bu sonuç doğrudan karşılıklılıktaki simülasyonun aksine verilen eşitsizliğe uygun çıktı. \(q \) değeri yükseldikçe işbirliği daha etkili bir strateji olacak mı diye bir uç değeri de denedim.

\(q = 0.99 \)

Gerçekten de bu simülasyona göre ihanet stratejisi uygulayanlar popülasyonda tanındığı zaman toplum işbirliği stratejisini benimsemeye başlıyor.

İki kişilik oyun tipleri

Evrimsel oyun teorisi konusunda okuduğum son kitapta nedense fitness generating function kısmını bir türlü anlayamadım. Teknik kitaplar her şeyi oldukça basit bir şekilde anlatmalarıyla meşhur değiller zaten. Bunun üzerine kitapta bahsedilen makalelere bakmaya başladım ve bu da bir yerde çıkmaz sokakla bitti. Akademisyen olmadığım için belli bir altyapımın olmadığının farkındayım ama bu altyapıyı sağlamanın kolay bir yolu da olmalıdır diye düşünüyordum. Neyse, aramaya devam edeceğim. Bu sırada yoluma Jun Tanimoto’nun Fundamentals of Evolutionary Game Theory and its Applications adlı kitabından devam etmeye karar verdim.

Kitap önce iki kişilik oyunları belli sınıflara ayırıp biraz teoriden bahsediyor. Benim hedefim kuru teori yerine bu oyunları, algoritmaları ve simülasyonlarını programlamak olduğundan hemen python denemelerine başladım.

Bu denemelerden bahsetmeden önce biraz notasyondan bahsetmem gerekecek. İki kişilik ve iki stratejili oyunlarla başladım. Buna 2×2 oyunlar da deniyor, payoff matrisimiz de 2×2 matris.

OyuncularCooperation (C)Defection (D)
Cooperation (C)R, RS, T
Defection (D)T, SP, P

Kısaca bu tabloyu kullanmayı anlatayım. Cooperation işbirliği stratejisi oluyor, Defection da ihanet. Her adımda seçilen iki kişi birbiriyle stratejilerine göre bir oyun oynuyor. Bütün oyun bu tabloda kodlanmış durumda ve seçilen stratejilere göre kazançları da R, S, T, P değerleriyle tanımlanıyor. Oyuncuların biri soldaki sütundaki stratejilerden birine sahip, diğeri de ilk satırdaki stratejilerden birini oynuyor. Bu durumda şu eşleşmeler mümkün:

Birinci oyuncu: Cooperation   İkinci oyuncu: Cooperation
Birinci oyuncu: Cooperation   İkinci oyuncu: Defection
Birinci oyuncu: Defection     İkinci oyuncu: Cooperation
Birinci oyuncu: Defection     İkinci oyuncu: Defection 

Her eşleşme için tablodaki eşleşme hücresindeki değerler oyuncuların kazançlarını belirliyor. Bu durumda yukarıdaki eşleşmeler için sonuçlar şöyle olacak:

Birinci oyuncu: Cooperation   İkinci oyuncu: Cooperation
Bu durumda kazançlar tablodaki R, R değerlerinin olduğu hücreden okunacak. Yani birinci oyuncunun kazancı R değeri kadar, ikinci oyuncunun kazancı da R değeri kadar olacak.

Birinci oyuncu: Cooperation   İkinci oyuncu: Defection
Bu durumda kazançlar tablodaki S, T değerlerinin olduğu hücreden okunacak. Yani birinci oyuncunun kazancı S değeri kadar, ikinci oyuncunun kazancı da T değeri kadar olacak.

Birinci oyuncu: Defection     İkinci oyuncu: Cooperation
Bu durumda kazançlar tablodaki T, S değerlerinin olduğu hücreden okunacak. Yani birinci oyuncunun kazancı T değeri kadar, ikinci oyuncunun kazancı da S değeri kadar olacak.

Birinci oyuncu: Defection     İkinci oyuncu: Defection
Bu durumda kazançlar tablodaki P, P değerlerinin olduğu hücreden okunacak. Yani birinci oyuncunun kazancı P değeri kadar, ikinci oyuncunun kazancı da P değeri kadar olacak.

Programda kullanırken bu tabloyu kitapta da gösterildiği gibi aşağıdaki M matrisi haline getirdim.

\(M = \begin{bmatrix} R&&S\\T&&P \end{bmatrix} \)

Kitapta olası oyunları incelemek için şu iki değer de hesaplanıyor.

\(D_{g} = T - R \)
\(D_{r} = P - S \)

\(D_{g} \) terimi sıfırdan büyükse davranış risk alma ikilemine sahip oluyormuş. \(D_{r} \) terimi sıfırdan büyükse bu sefer de davranış riskten kaçma ikilemine sahip oluyormuş.

Önce bu sınıflamayı kısaca vereyim, ilerki yazılarda bu sınıflamalar üzerine biraz daha düşünmeye çalışırım.

Tutsak ikilemi: Burada hem \(D_{g} \) hem de \(D_{r} \) sıfırdan büyük

Tavuk: Bu tip oyunlarda \(D_{g} \) sıfırdan büyük ve \(D_{r} \) sıfırdan büyük değil.

Geyik avı: Bu tip oyunlarda sadece \(D_{r} \) sıfırdan büyük oluyor.

İkilemsiz tip: Bu tip oyunlarda ne \(D_{g} \) ne de \(D_{r} \) sıfırdan büyüktür.

Bu oyun türleri hakkında internette bir sürü bilgi bulunabilir. Bu nedenle kitapta verilen algoritma ve bu tipler için yazdığım programı göstereceğim. Son olarak da bu tiplerin sonuçlarını grafik üzerinden kısaca anlatmaya çalışacağım.

Seçilen oyuncunun stratejisi iki ihtimalden (Cooperation ya da Defection) biri olacak. Bunları da şu iki vektörle tanımlayalım:

\(\vec{e_{1}} = \begin{pmatrix}1\\0\end{pmatrix} \) Cooperation stratejisini oynayan bireyi simgelesin.

\(\vec{e_{2}} = \begin{pmatrix}0\\1\end{pmatrix} \) de Defection stratejisini oynayan bireyi simgelesin.

Stratejilerin nüfustaki oranları da aşağıdaki vektörle (Transpose edilmiş hali) verilsin.

\(s^{T} =\) $latex \begin{pmatrix}s_1&&s_2\end{pmatrix}

Cooperation stratejisini kullanan bir bireyin rastgele seçilen başka bir bireye karşı beklenen kazancı şu şekilde veriliyor:

\(e_{1}^{T} \cdot M \cdot s \)

Aynı şekilde Defection stratejisini kullanan bir bireyin rastgele seçilen başka bir bireye karşı beklenen kazancı da şu şekilde olur:

\(e_{2}^{T} \cdot M \cdot s \)

Bireylerin çoğalma mekanizması da şu şekilde tanımlanıyor:

\(\frac{\overset {.}{s_{i}}}{s_{i}} = {e_{i}}^{T} \cdot M \cdot s – {s}^{T} \cdot M \cdot s \)

Burada \(\overset {.}{s_{1} \) terimi Cooperation stratejisini kullanan bireylerin oranının değişim miktarını gösteriyor. \(\overset {.}{s_{2} \) terimi de Defection stratejisini kullanan bireylerin oranının değişim miktarını gösteriyor. \(s_{1} \) ve \(s_{2} \) terimleri de bu iki stratejiyi kullanan bireylerin toplam nüfustaki oranları ve toplamları da haliyle bir olmak zorunda. Bu iki toplam hep sabit kaldığından birinin oranındaki artış ile diğer grubun oranındaki azalış da büyüklük olarak birbirine eşit olmak zorunda.

Sistemin dinamik davranışını incelemek için de aşağıdaki python programını yazdım. Program açıklamalarını kodun içine eklemeye çalıştım. İstenen oyun tipi simülasyonunu elde etmek için payoff değişkeni program başındaki tiplerden birine eşitlemek yeterli olacaktır.

import numpy as np
import matplotlib.pyplot as plt

strategies = np.array([[1, 0], [0, 1]]) #strateji vektörleri

prisoner_dilemma_payoff = np.array([[5, 1], [7, 3]]) #tutsak ikilemi oyun tipi

chicken_dilemma_payoff = np.array([[5, 1], [7, 0]]) #tavuk oyun tipi

stag_hunt_dilemma_payoff = np.array([[8, 1], [7, 3]]) #geyik avı oyun tipi

trivial_payoff = np.array([[7, 3], [5, 1]]) #ikilemsiz oyun tipi

number_of_iterations = 100 # her bir oyunun kaç jenerasyon boyunca oynanacağı

payoff = prisoner_dilemma_payoff  # seçilen oyun tipi

increment = 0.01
steady_states = np.array([[0, 0]]) # oyun bittiğinde nüfusun hangi oranlardan oluştuğu

for s1 in np.arange(0, 1.0, increment): #0-1 aralığında değişik başlangıç  konumları yaratılıyor
    s2 = 1 - s1  # nüfus oranları toplamı 1 olmalı
    species = np.array([s1, s2]) # ilk durum olarak türlerin oranı belirlenmiş oldu

    for i in range(0, number_of_iterations): # bu şartlar için oyun belirlenen jenerasyon miktarı kadar oynanacak

        difference = (strategies.dot(payoff).dot(species) - species.dot(payoff).dot(species))*species # strateji vektörlerinin değişim miktarları hesaplanıyor
        species = species + difference # bu stratejileri kullanan bireylerin nüfustaki oranları güncelleniyor
        species = np.clip(species, 0, 1) # nüfus oranlarının her zaman 0 ile 1 arasında olması sağlanıyor

    steady_states = np.append(steady_states, np.array([[s1, species[0]]]), axis = 0)  # oyunun sonunda elde edilen oranlar kayıtlara ekleniyor

#Alttaki kısımda da simülasyon sonuçları grafik olarak hazırlanıyor
plt.plot(steady_states[:, 0], steady_states[:, 1])
plt.ylabel("1. türün sondaki oranı")
plt.xlabel("1. türün başlangıçtaki oranı")
plt.show()
Tutsak ikilemi oyununda başlangıçta Cooperation stratejisini uygulayanların oranı ne olursa olsun oyunun sonunda nüfus tamamen Defection stratejili bireyler kalıyor.

Tutsak ikileminde kullandığım matris \(M = \begin{bmatrix} 5&&1\\7&&3 \end{bmatrix} \) şeklindeydi.

Tavuk tipi oyunda dahemen hemen her nüfus dağılımı için iki stratejinin de temsil edildiği denge durumları meydana geldi.

Tavuk tipi oyunda kullandığım matris de \(M = \begin{bmatrix} 5&&1\\7&&0 \end{bmatrix} \) şeklindeydi.

Geyik avı tipi oyunda bir orana kadar hep Cooperation stratejisi kazanıyor, o orandan sonra ise Defection. İki stratejinin bir arada süregeldiği bir durum olmuyor.

Geyik avında kullandığım matris de \(M = \begin{bmatrix} 8&&1\\7&&3 \end{bmatrix} \) şeklindeydi.

İkilemsiz oyun tipinde de her dağılım için Cooperation stratejisi nüfusu ele geçiriyor.

Bu oyun tipinde de \(M = \begin{bmatrix} 7&&3\\5&&1 \end{bmatrix} \) matrisini kullandım.

Şimdilik bu yazıyı burada bitireyim. Amacım basit oyunların simülasyonlarını kolayca yapmanın bir yolunu bulmaktı ve şimdilik bu kitaptaki yöntemle çalışabiliyorum gibi gözüküyor. Umarım konular ilerledikçe buna yakın bir performans gösterebilirim.

Tüketici-Kaynak modeli

Bu model kitapta çok kısa anlatılmıştı. Bu nedenle kaynak olarak verilen makaleleri bulup anlamaya çalıştım. Bu modelde de iki tür kullandım. Bunlar iki bitki türü ve kaynak olarak da topraktaki besleyici maddeler ve güneş ışığı var.

Önce türlerin gelişimini kontrol eden denklemleri vereyim. İtiraf edeyim, bu denklemler önceki modellerdeki denklemlerden daha karmaşık. İki türün seçtiği stratejiler aşağıdaki \(A_1\) ve \(A_2\) parametreleriyle seçiliyor.

\(\frac{dB_1}{dt} = B_1 \cdot MIN(\frac{r \cdot N \cdot A_1}{N + k_N} – R – d , \frac{r \cdot L (1 – A_1)}{L + k_L} – R – d) \)

\(\frac{dB_2}{dt} = B_2 \cdot MIN(\frac{r \cdot N \cdot A_2}{N + k_N} – R – d , \frac{r \cdot L (1 – A_2)}{L + k_L} – R – d) \)

Bu denklemlerde MIN işlemi, iki değerin en küçüğünü seçiyor. Görüldüğü gibi modelde bir sürü sabit ve parametre var. Bunları kısaca açıklamaya çalışayım.

\(B_1\): Birinci türün toplam biyokütlesi. Burada bitki türlerinden bahsettiğimiz için nüfus yerine toplam ne kadar biyokütle oluşturduğuna bakılmış denklemlerde. Her iterasyonda bu değer yeniden hesaplanacak.

\(B_2\): İkinci türün toplam biyokütlesi. Bu değer de her iterasyonda yeniden hesaplanıyor.

Şimdi bu biyokütleleri hesaplamak için kullanacağım diğer parametreleri açıklayayım.

\(r \): Her bitkinin maksimum büyüme oranı. Bitkilerin biyokütlelerinin artışı bu sabitle doğru orantılı, yani bu sabit ne kadar büyükse bitkilerin büyüme oranı da o kadar hızlı. Doğru orantıyı görmek için ilk denklemde önce \(B_1 \cdot {r \cdot N \cdot A_1}{N + k_N} \) kısmına bakalım. Burada \(r \) sayısı büyürse çarpım da büyür. Aynı şekilde minimum işleminin diğer adayında da \(r \) sayısı pay kısmında bulunuyor. \(B_1 \cdot r \cdot L (1 – A_1) \). Yani \(r \) sayısı büyüdükçe eşitliğin sağ tarafı da büyüme eğiliminde.

\(A_1\): Birinci tür için köke ayrılan biyokütle. Bu parametre birinci türün seçtiği stratejiyi de belirliyor. Yani daha çok besine mi güneş ışığına mı ağırlık veriyor.

\(A_2\): İkinci tür için köke ayrılan biyokütle. Bu parametre ikinci türün seçtiği stratejiyi de belirliyor. Yani daha çok besine mi güneş ışığına mı ağırlık veriyor.

\(k_N\): Topraktaki besine dayalı büyüme için yarı doyma sabiti. Bu sabit paydada olduğu için arttıkça topraktaki besinden gelen biyokütle artışı azalıyor. Bu sabit sayesinde topraktaki besin değişikliğinin ani etkileri biraz azaltılmış oluyor. Yani \(N \) değeri aniden büyürse ya da küçülürse türün biyokütlesi daha az bir oranda büyüyüp küçülecek. Bunu görmek için \(\frac{N}{N + k_N} \) ifadesiyle biraz oynamak yeterli.

\(k_L\): Güneş ışığına bağlı büyüme için yarı doyma sabiti. Bu sabit paydada olduğu için arttıkça güneş ışığından gelen biyokütle artışı azalıyor. Bu sabit sayesinde güneş ışığındaki değişikliğinin ani etkileri biraz azaltılmış oluyor. Yani \(L \) değeri aniden büyürse ya da küçülürse türün biyokütlesi daha az bir oranda büyüyüp küçülecek. Bunu görmek için \(\frac{L}{L + k_L} \) ifadesiyle biraz oynamak yeterli.

\(R \): Solunum oranı. Türlerden bağımsız bir sabit. Sanırım solunum sırasında besinlerin yakılması nedeniyle biyokütledeki bir tür azalmayı modelliyor. Bu sabitle biyokütlenin çarpımı çıkarma işlemi yüzünden biyokütleyi her zaman azaltma eğiliminde.

\(d \): Kayıp oranı. Türlerden bağımsız ve iç ya da dış nedenlerle herhangi bir şekilde ölümlerin modellendiği bir sabit. Bu sabitle biyokütlenin çarpımı o anki biyokütleden çıkarıldığı için beklendiği gibi biyokütleyi her zaman azaltma eğiliminde.

\(N \): Topraktaki toplam besin miktarı. Bu modelde kullanılan kaynaklardan birisi bu. Bitkiler çoğaldıkça topraktaki besin miktarı azalacaktır ve bitkiler ölünce toprak da besin bakımından zenginleşecektir. Bu kaynağın zamana göre değişim denklemini aşağıda vereceğim.

\(L \): Bitkiler tarafından kullanılabilen güneş ışığı miktarı. Bu da bitkilerin kullanabildiği kaynaklardan ikincisi.

Şimdi kaynakların değişimini modelleyen denklemlere bakalım.

\(\frac{dN}{dt}=a (T – N – B_1 \cdot p – B_2 \cdot p) \)

\(– MIN(\frac{r\cdot N \cdot A_1}{N + k_N} – R, \frac{r\cdot L (1 – A_1}{L + k_L} – R) \cdot B_1 \cdot p \)

\(– MIN(\frac{r\cdot N \cdot A_2}{N + k_N} – R, \frac{r\cdot L (1 – A_2}{L + k_L} – R) \cdot B_2 \cdot p \)

\(L = \frac {L_0}{1 + \alpha \cdot B_1(1-A_1) + \alpha \cdot B_2(1 – A_2)} \)

Bu iki denklemin parametrelerini açıklayayım:

\(T \) : Habitattaki toplam besin miktarı. Bu toplama topraktaki, ölü bitkilerdeki ve bitki biyokütlelerindeki besinlerin hepsi dahil.

\(p \): Bitki dokularındaki besin konsantrasyonu.

\(a \): Mineralleşme oranı. Yani ölü biyokütlenin hangi oranda besine dönüştüğünü belirliyor. \(a (T – N – B_1 \cdot p – B_2 \cdot p) \) ifadesindeki \(T – N – B_1 \cdot p – B_2 \cdot p \) kısmına bakarsak önce habitattaki toplam besinden topraktaki besin miktarını çıkarıyoruz. Sonra kalandan birinci türün biyokütlesiyle bitki dokusundaki besin konsantrasyonunu yani birinci türdeki toplam besin miktarını çıkarıyoruz. Sonra yine kalandan ikinci türdeki toplam besin miktarını çıkarıyoruz. Dolayısıyla elimizde kalan miktar aslında ölü bitkilerdeki besin miktarını verir. Bu miktarıda mineralleşmeyle çarparsak o iterasyonda ölü bitkilerden toprağa ne kadar besin geri dönüşü olacağını buluruz. Tabii ki mineralleşme oranının kullanıldığı ilk kaynak denkleminde minimum işleminin yapıldığı bir terim daha var. Dikkat edersek bu terim neredeyse biyokütle değişim denklemlerinin aynısı. Farkların biri kayıp parametresi burada yok, diğeri de biyokütledeki artışın bitki besin konsantrasyonuyla çarpılması. Bu iki fark aslında bu terimin biyokütledeki değişim sayesinde ne kadar besinin topraktan uzaklaştırıldığını gösteriyor ve bu nedenle topraktaki besin miktarına etkisi de çıkarma işleminden görüldüğü gibi negatif.

\(L_0\): Güneş ışığı sabiti.

\(\alpha \): Birim yaprak biyokütlesi başına soğrulan ışık miktarı. Bu sabit paydada olduğundan büyüdükçe ortamdaki diğer bitkilere kalan ışık miktarı da azalıyor. Paydadaki diğer terimler de toplam ışıktan faydalanan biyokütle olduğunu görmek kolay. Yani bütün ifade bir tür biyokütle başına düşen ışık miktarı gibi düşünülebilir.

Bu modeli denemek için aşağıdaki basit kodu yazdım. Sabitlerle ve stratejilerle oynayarak çok değişik çoğalma davranışları görmek mümkün. Değişkenleri yazıdaki şekliyle tanımladım ki, değişiklik yapmak daha kolay olsun.

import numpy as np
import matplotlib.pyplot as plt

B1 = 0.5  # Biomass of population 1
B2 = 0.5  # Biomass of population 2
N = 1.0  # available soil nutrient
L = 1.0  # light availability
A1 = 0.8  # fraction of biomass allocated to root by population 1
A2 = 0.2  # fraction of biomass allocated to root by population 2
r = 5.0  # per capita maximal rate of plant growth
kn = 0.1  # 1/2 saturation constant for nutrient
kl = 0.1  # 1/2 saturation constants for light
R = 0.2  # density independent per capita respiration rate
d = 0.1  # density independent per capita loss rate


def G1(L, N):
    return r*N*A1/(N+kn) - R - d


def G2(L, N):
    return r*L*(1-A2)/(L+kl) - R - d


T = 3.0  # Total soil nutrient in habitat
a = 0.5  # mineralization rate
p = 0.4  # plant tissue nutrient concentration

L0 = 2.0  # solar constant
alpha = 0.01  # light decay rate per unit leaf biomass

populations = np.empty((0, 4), int)
populations = np.append(populations, np.array(
    [[B1, B2, N, L]]), axis=0)

for i in range(1, 100):
    B1 = B1 + B1*min(G1(L, N), G2(L, N))
    if B1 < 0:
        B1 = 0
    B2 = B2 + B2*min(G1(L, N), G2(L, N))
    if B2 < 0:
        B2 = 0

    N = N + a*(T - N - B1*p - B2*p) - min(G1(L, N), G2(L, N)) * \
        B1*p - min(G1(L, N), G2(L, N))*B2*p
    if N < 0:
        N = 0
    L = L0 / (1 + alpha*B1*(1 - A1) + alpha*B2*(1 - A2))
    if L < 0:
        L = 0

    populations = np.append(populations, np.array(
        [[B1, B2, N, L]]), axis=0)


f, (ax1, ax2, ax3, ax4) = plt.subplots(4)
line1, = ax1.plot(populations[:, 0], color="b")
line2, = ax2.plot(populations[:, 1], color="r")
line3, = ax3.plot(populations[:, 2], color="g")
line4, = ax4.plot(populations[:, 3], color="y")

ax1.set_ylabel("B1")
ax2.set_ylabel("B2")
ax3.set_ylabel("Besin")
ax4.set_ylabel("Güneş")
ax2.set_xlabel("zaman")
plt.show()

Bu da yukarıdaki programın çıktısı

Lotka-Volterra modeli

Oyun teorisi her zaman ilgimi çekmiştir. Geçtiğimiz yıllarda evrimsel oyun teorisine bir giriş yaptım fakat kullandığım kitapta birçok şey başlangıç seviyesi için uygun olmadığından çok ileri gidemedim. Geçenlerde başka bir kitap buldum ve şansımı bir de bununla deneyeyim dedim.

İlk pratik konular popülasyon modelleme üzerineydi. Bu modellerden birini bu yazıda anlamaya ve anlatmaya çalışacağım.

Bu modelde iki türden oluşan popülasyonların bir türü ele alınıyor. En basit örnekleri av ve avcılardan şeklindeki türler. Avcıların doğum oranı çevredeki av oranına bağlı. Ne kadar çok avlanabilirlerse o kadar iyi beslenirler ve üremeleri kolaylaşır. Avcıların ölüm miktarı ise normal şartlarda sadece kendi nüfuslarına bağlı. Aynı şekilde avların üreme miktarı sadece kendi nüfuslarına bağlıyken ölüm sayıları tabii ki avcıların sayısına da bağlı. Ortamdaki kaynakların hem av hem de avcı için yeterli olduğunu varsayalım.

Bu düşünceler doğrultusunda şimdi Lotka-Volterra denklemlerini yazayım. x av türünün nüfusu, y de avcı türünün nüfusu olsun.

\(\frac{dx}{dt} = \alpha x – \beta x \cdot y \)

\(\frac{dy}{dt} = \delta x \cdot y – \gamma y \)

Bu modeli biraz açıklamaya çalışayım. \(\frac{dx}{td} \) terimi av türünün nüfusunun zamana göre değişimini göstermekte. Yani bir süre geçince bu nüfusun artışını ya da azalışını, kısaca değişimini gösteriyoruz. Bu değişim o anki nüfus sayılarıyla ve bası sabitlerle orantılı. Aynı şekilde \(\frac{dy}{td} \) terimi de avcı nüfusundaki zamana göre değişimi ifade etmekte.

Denklemlerdeki sabitlere bakalım şimdi: Öncelikle bütün sabitlerin pozitif sayılar olduğunu varsayalım. Bu sayede o sabitlere daha alışık olduğumuz anlamlar yükleyebileceğiz.

\(alpha \) : av türünün doğum oranı. Doğum oranı av türünde sadece av türünün popülasyonunu etkiliyor. Bunu \(\alpha x \) ifadesinden görebiliyoruz. Av türünün çoğalması için çevredeki kaynaklara ve kendi türüne ihtiyacı var ve çevredeki kaynakların yeterli olduğunu varsaymıştık.

\(\beta \): bunu av türündeki ölüm oranı diye yorumlayabiliriz. Genelde iki tür arasındaki etkileşim faktörü olarak kullanılıyor ama bu denklemde ölüm oranı olarak çalıştığını görebiliriz. Öncelikle $\beta x \cdot y $ terimi diğer terimden çıkarılıyor, yani bir azalmayı simgeliyor. Diğer taraftan iki popülasyonun sayısına da etki ediyor. Yani bu sayılar arttıkça etkisi de artıyor. Bunu avcı sayısı arttıkça bu avcıların avlayacağı av sayısının artması şeklinde yorumlayabiliriz. Aynı şekilde av sayısı artarsa da avcılar daha kolay av bulacağından, az avcı olsa bile yine de daha çok av ölecektir.

\(\delta \): Bu sabiti de avcıların doğum oranı olarak düşünebiliriz. Aslında \(\beta \) sabiti gibi iki tür arasındaki bir etkileşim sabiti ama \(\delta x \cdot y \) teriminin pozitif olduğuna bakarsak avcı popülasyonunu artıracak bir etki yaptığını görürüz. Ayrıca iki popülasyona da bağlı bir artış sağlıyor. Eğer avcı sayısı artarsa türün üremesi kolaylaşacak. Eğer av sayısı artarsa avcılar daha iyi beslenebilecek ve yine daha çok üreyebilecekler.

\(\gamma \): Bu sabiti de avcıların ölüm oranı şeklinde düşünebiliriz. \(\gamma y \) teriminin diğer terimden çıkarılması popülasyonda azaltıcı bir etki yapacaktır ve bu da ölümü çağrıştırabilir. Av türü avcılar için tehdit olmadığında bu sabit av nüfusundan bağımsız bir etki yapmakta.

Şimdi bu denklemleri python programlarıyla deneyerek popülasyonların nasıl değiştiğinin grafiklerini elde etmeye çalışacağım.

import numpy as np
import matplotlib.pyplot as plt

x = 10
y = 10
alpha = 0.1
beta = 0.04
gamma = 0.04
delta = 0.01

populations = np.empty((0, 2), int)
populations = np.append(populations, np.array(
    [[x, y]]), axis=0)

for i in range(1, 1000):
    x = x + alpha*x - beta*x*y
    if x < 0:
        x = 0
    y = y + delta*x*y - gamma*y
    if y < 0:
        y = 0
    populations = np.append(populations, np.array(
        [[x, y]]), axis=0)

f, (ax1, ax2) = plt.subplots(2)

line1, = ax1.plot(populations[:, 0], color="b")
line2, = ax2.plot(populations[:, 1], color="r")

ax1.set_ylabel("Av")
ax2.set_ylabel("Avcı")
ax2.set_xlabel("zaman")
plt.show()

Av sayısı artışından hemen sonra bol besin bulan avcıların sayısı da hızla yükselmiş. Avcı sayısının yükselmesiyle av nüfusunda aşırı avlanma nedeniyle hızlı bir düşüş görülüyor. Grafikten ayrıca bu davranışların sürekli tekrarlandığını görüyoruz. Bu grafiğe bakıp her zaman bu tür bir nüfus değişimi gözlemleneceği düşünülmemeli ama. Sabit sayılara ve ilk durumdaki nüfus değerlerine göre çok farklı bir davranış da görülebilir. Şimdi buna bir örnek vereyim.


import numpy as np
import matplotlib.pyplot as plt

x = 10
y = 10
alpha = 0.1
beta = 0.1092
gamma = 0.04
delta = 0.01

populations = np.empty((0, 2), int)
populations = np.append(populations, np.array(
    [[x, y]]), axis=0)

for i in range(1, 5000):
    x = x + alpha*x - beta*x*y
    if x < 0:
        x = 0
    y = y + delta*x*y - gamma*y
    if y < 0:
        y = 0
    populations = np.append(populations, np.array(
        [[x, y]]), axis=0)

f, (ax1, ax2) = plt.subplots(2)

av, = ax1.plot(populations[:, 0], color="b")
avci, = ax2.plot(populations[:, 1], color="r")

ax1.set_ylabel("Av")
ax2.set_ylabel("Avcı")
ax2.set_xlabel("zaman")
plt.show()

Bu programda sadece \(\beta \) sabitini yükselttim, yani avcıların avlanma oranını artırdım. Bir noktadan sonra aşırı avlanmadan sonra av türünün nüfusu tamamen tükeniyor. Bu noktadan sonra avcılar da yok olmaya başlıyor.

Zamanı bu grafikte daha da ileriye götürdüm ki türlerin bir daha dirilmediği daha iyi görülebilsin. Demek ki seçilen sabitlere göre iki türün de ortadan kalkması mümkün.

Peki av türü ilk azalmadan sonra nasıl yeniden çoğalmayı başardı? Bunu görebilmek için grafiği çok büyütmek lazım.

Bu grafikte sadece av nüfusunu 3. ile 9. iterasyon arasında büyüttüm. Avcı grafiğine dokunmadım. Burada görüyoruz ki av nüfusu sıfıra çok yaklaşmış (\(5\cdot 10^{-11}\) ). Evet bu sayı gerçek hayatta çok anlamsız ama modeli tanımlayan diferansiyel denklemler için geçerli ve yeterli sayılar.

Bir de bu iki türün popülasyonlarını zamana göre değil de birbirlerine göre gösterebiliriz. Buna faz uzayı deniyor sanırım, ingilizcesi phase-space plot. Programı bir de bu gösterim için düzenledim.


from random import betavariate
import numpy as np
import matplotlib.pyplot as plt

x = 10
y = 10
alpha = 0.1
beta = 0.04
gamma = 0.04
delta = 0.01

populations = np.empty((0, 2), int)
populations = np.append(populations, np.array(
    [[x, y]]), axis=0)

for i in range(1, 1000):
    x = x + alpha*x - beta*x*y
    if x < 0:
        x = 0
    y = y + delta*x*y - gamma*y
    if y < 0:
        y = 0
    populations = np.append(populations, np.array(
        [[x, y]]), axis=0)

plt.plot(populations[:, 0], populations[:, 1])
plt.ylabel("avcı")
plt.xlabel("av")
plt.show()

Bu da programı çalıştırınca çıkan grafik:

Bu model çevremizdeki türleri çok iyi modellemese de uzun süre ekonomi alanında kullanılmış.