先引用一下果殼網對謝林模型的介紹

該模型由馬里蘭大學的經濟學家托馬斯·謝林建立。主要感興趣的隔離有兩種：一種是種族隔離，另一種是由收入差異導致的隔離。

• 種族隔離：

每一個紅點代表城市中白人的聚集地，每一個藍點代表城市中美國黑人的聚集地，每一個黃點代表城市中拉美人的聚集地，綠點則代表城市中亞裔的聚集地。從這幅圖可以看到非常明顯的種族隔離。

• 收入隔離：

圖片每個紅點代表富人，淺藍色的點代表窮人，而深藍點則代表中產階層。從這幅圖可以看到收入隔離現象儘管非常明顯，但並沒有種族隔離現象明顯。

謝林構建了一個基於代理人的模型。該類模型要點有三個：代理人，規則，以及群體表現。

謝林模型是個關於人們選擇在哪居住的模型。他把每個人放置在棋盤上，整個城市看作是一塊巨大的棋盤。棋盤上每一個小格子可以住一個人也可以空著。

假設紅格子表示住著富人，灰格子表示住了窮人，白格子表示沒有主人。對在中間X格子的富人來說，周圍7個鄰居有3個是同類，他應該留下還是搬走？

謝林設定了規則稱作閾值基本準則：每個人都有一個閾值（threshold），根據閾值做決定是留下還是搬走。比如假設X的閾值是33%，則此時他會選擇留下。而當他的富人鄰居有一位搬走時，由於只剩2/7是同類所以那時會選擇搬走。

居住地總數：60*60；居民總數：2800

滿意度閾值：(0,1) 設定

紅藍兩色代表不同膚色，綠色代表空

每個節點相似比低於相似度閾值隨機搬家，所有節點相似比不低於相似度閾值，停止迭代

下面是我用python實現的模擬：python2.7+matplotlib+numpy

import random
import numpy as np
import matplotlib.pyplot as plt
import time
# plt.ion()
vector = range(1,3601) #1-3600
random.shuffle(vector)
white = vector[0:1400]
black = vector[1400:2800]
temp = vector[2800:3600]
tmap = np.zeros((62,62,3))
# init
for i in range(1,3601):
    if i in white:
        tmap[(i - 1) / 60 + 1][(i - 1) % 60 + 1] = [1, 0, 0]
    elif i in black:
        tmap[(i - 1) / 60 + 1][(i - 1) % 60 + 1] = [0, 0, 1]
    else:
        tmap[(i - 1) / 60 + 1][(i - 1) % 60 + 1] = [0, 1, 0]
row = [-1,-1,-1, 0, 0, 1, 1, 1] #neiborhood
col = [-1, 0, 1,-1, 1,-1, 0, 1]
simv = []
for m in range(7):
    thd = 0.3+m*0.1  #threshold
    print 'Threshold:'+str(thd)
    sim = 0.0 #similarity
    nosat = 1
    times = 0
    while nosat and times<50:
        times=times+1
        sim = 0.0
        bnosat = 0
        wnosat = 0
        nosatv=[]
        for i in range(1,61):
            for j in range(1,61):
                if(tmap[i][j][1] != 1):
                    ns=0.0
                    n=0.0
                    for k in range(8):
                        if tmap[i+row[k]][j+col[k]][0]==1 or tmap[i+row[k]][j+col[k]][2]==1:
                            n=n+1.0
                            if tmap[i+row[k]][j+col[k]][0] == tmap[i][j][0]:
                                ns=ns+1.0
                    if n!=0:
                        sim = sim + ns/n #similarity plus
                        
                        if(ns/n < thd):
                            if(tmap[i][j][0]==1):
                                wnosat = wnosat+1
                            else:
                                bnosat = bnosat+1
                            nosatv.append([i, j])
                else:
                    nosatv.append([i, j])
        nosat = wnosat + bnosat
        print 'Unsatisfied:' + str(nosat / 2800.0*100)+'%'
        #-------------move--------------
        random.shuffle(nosatv)
        wnosatv=nosatv[0:wnosat]
        bnosatv=nosatv[wnosat:wnosat+bnosat]
        tempv = nosatv[wnosat+bnosat:]
        for i in range(wnosat):
            tmap[wnosatv[i][0]][wnosatv[i][1]]=[1,0,0]
        for i in range(bnosat):
            tmap[bnosatv[i][0]][bnosatv[i][1]]=[0,0,1]
        for i in range(tempv.__len__()):
            tmap[tempv[i][0]][tempv[i][1]] = [0, 1, 0]

    sim = sim/2800.0
    simv.append(sim)
    print 'Similarity:' + str(sim*100)+'%'

    plt.figure()
    plt.axis("off")
    plt.imshow(tmap, interpolation="nearest")
    plt.title('Threshold = '+str(thd)+'    Similarity:'+ str(sim*100)+'%')
    plt.show()

plt.figure()
plt.plot([0.3+m*0.1 for m in range(7)],simv,'-')
plt.title('Similarity of different threshold')
plt.xlabel('threshold')# make axis labels
plt.ylabel('Similarity')
plt.xlim(0.2, 1.0)
plt.xlim(0.2, 1.0)
plt.show()