【元胞自动机】基于matlab元胞自动机模拟SEIR传播模型【含Matlab源码 2156期】
一、元胞自動機(jī)簡介
1 元胞自動機(jī)發(fā)展歷程
最初的元胞自動機(jī)是由馮 · 諾依曼在 1950 年代為模擬生物 細(xì)胞的自我復(fù)制而提出的. 但是并未受到學(xué)術(shù)界重視.1970 年, 劍橋大學(xué)的約翰 · 何頓 · 康威設(shè)計了一個電腦游戲 “生命游戲” 后, 元胞自動機(jī)才吸引了科學(xué)家們的注意.1983 年 S.Wolfram 發(fā)表了一系列論文. 對初等元胞機(jī) 256 種 規(guī)則所產(chǎn)生的模型進(jìn)行了深入研究, 并用熵來描述其演化行 為, 將細(xì)胞自動機(jī)分為平穩(wěn)型, 周期型, 混沌型和復(fù)雜型.
2 對元胞自動機(jī)的初步認(rèn)識
元胞自動機(jī)(CA)是一種用來仿真局部規(guī)則和局部聯(lián)系的方法。典型的元胞自動機(jī)是定義在網(wǎng)格上的,每一個點(diǎn)上的網(wǎng)格代表一個元胞與一種有限的狀態(tài)。變化規(guī)則適用于每一個元胞并且同時進(jìn)行。典型的變化規(guī)則,決定于元胞的狀態(tài),以及其( 4 或 8 )鄰居的狀態(tài)。
3 元胞的變化規(guī)則&元胞狀態(tài)
典型的變化規(guī)則,決定于元胞的狀態(tài),以及其( 4 或 8 )鄰居的狀態(tài)。
4 元胞自動機(jī)的應(yīng)用
元胞自動機(jī)已被應(yīng)用于物理模擬,生物模擬等領(lǐng)域。
5 元胞自動機(jī)的matlab編程
結(jié)合以上,我們可以理解元胞自動機(jī)仿真需要理解三點(diǎn)。一是元胞,在matlab中可以理解為矩陣中的一點(diǎn)或多點(diǎn)組成的方形塊,一般我們用矩陣中的一點(diǎn)代表一個元胞。二是變化規(guī)則,元胞的變化規(guī)則決定元胞下一刻的狀態(tài)。三是元胞的狀態(tài),元胞的狀態(tài)是自定義的,通常是對立的狀態(tài),比如生物的存活狀態(tài)或死亡狀態(tài),紅燈或綠燈,該點(diǎn)有障礙物或者沒有障礙物等等。
6 一維元胞自動機(jī)——交通規(guī)則
定義:
6.1 元胞分布于一維線性網(wǎng)格上.
6.2 元胞僅具有車和空兩種狀態(tài).
7 二維元胞自動機(jī)——生命游戲
定義:
7.1 元胞分布于二維方型網(wǎng)格上.
7.2 元胞僅具有生和死兩種狀態(tài).
元胞狀態(tài)由周圍八鄰居決定.
規(guī)則:
骷髏:死亡;笑臉:生存
周圍有三個笑臉,則中間變?yōu)樾δ?br /> 少于兩個笑臉或者多于三個,中間則變死亡。
8 什么是元胞自動機(jī)
離散的系統(tǒng): 元胞是定義在有限的時間和空間上的, 并且元 胞的狀態(tài)是有限.
動力學(xué)系統(tǒng): 元胞自動機(jī)的舉止行為具有動力學(xué)特征.
簡單與復(fù)雜: 元胞自動機(jī)用簡單規(guī)則控制相互作用的元胞 模擬復(fù)雜世界.
9 構(gòu)成要素
(1)元胞 (Cell)
元胞是元胞自動機(jī)基本單元:
狀態(tài): 每一個元胞都有記憶貯存狀態(tài)的功能.
離散: 簡單情況下, 元胞只有兩種可能狀態(tài); 較復(fù)雜情況下, 元胞具有多種狀態(tài).
更新: 元胞的狀態(tài)都安照動力規(guī)則不斷更新.
(2)網(wǎng)格 (Lattice)
不同維網(wǎng)格
常用二維網(wǎng)格
(3)鄰居 (Neighborhood)
(4)邊界 (Boundary)
反射型:以自己作為邊界的狀態(tài)
吸收型:不管邊界(車開到邊界就消失)
(5)規(guī)則(狀態(tài)轉(zhuǎn)移函數(shù))
定義:根據(jù)元胞當(dāng)前狀態(tài)及其鄰居狀況確定下一時刻該元胞狀態(tài)的動力學(xué)函數(shù), 簡單講, 就是一個狀態(tài)轉(zhuǎn)移函數(shù).
分類 :
總和型: 某元胞下時刻的狀態(tài)取決于且僅取決于它所有鄰居 的當(dāng)前狀態(tài)以及自身的當(dāng)前狀態(tài).
合法型: 總和型規(guī)則屬于合法型規(guī)則. 但如果把元胞自動機(jī) 的規(guī)則限制為總和型, 會使元胞自動機(jī)具有局限性.
(6)森林火災(zāi)
綠色:樹木;紅色:火;黑色:空地。
三種狀態(tài)循環(huán)轉(zhuǎn)化:
樹:周圍有火或者被閃電擊中就變成火。
空地:以概率p變?yōu)闃淠?br /> 理性分析:紅為火;灰為空地;綠是樹
元胞三種狀態(tài)的密度和為1
火轉(zhuǎn)化為空地的密度等于空地轉(zhuǎn)換為樹的密度(新長出來的樹等于燒沒的樹)
f是閃電的概率:遠(yuǎn)遠(yuǎn)小于樹生成的概率;T s m a x T_{smax}T smax
?是一大群樹被火燒的時間尺度
程序?qū)崿F(xiàn)
周期性邊界條件
購進(jìn)啊
其中的數(shù)字為編號
構(gòu)建鄰居矩陣
上面矩陣中的數(shù)字編號,對應(yīng)原矩陣相同位置編號的上鄰居編號,一 一對應(yīng)
同樣道理:
(7)交通概念
車距和密度
流量方程
守恒方程
時空軌跡(橫軸是空間縱軸為時間)
紅線橫線與藍(lán)色交點(diǎn)表示每個時間車的位置。
如果是豎線則表示車子在該位置對應(yīng)的時間
宏觀連續(xù)模型:
最常用的規(guī)則:
紅色條表示速度是滿的。
1 加速規(guī)則:不能超過v m a x ( 2 格 / s ) v_{max}(2格/s)v
max(2格/s)
2 防止碰撞:不能超過車距
理論分析:
結(jié)果分析: 密度與流量
第一個圖:橫坐標(biāo)是歸一化后的密度,縱坐標(biāo)是車流量。第二個圖:理論值與CA的結(jié)果
結(jié)果分析: 時空軌跡
中間的深色區(qū)域是交通堵塞的區(qū)域。
二、部分源代碼
% cellular automata
clc;clear;close all;
%%%%%% initialization of parameters
N = 100; % a half length of region
Dn = 0; % the number of death
rand(‘state’,121); % state generating random data
S = zeros(2N+1); % state matrix
r = rand(2N+1);
S(r>0.4) = 4; % uniform distribution of normal people before simulation
vc(1,:) = [1,1,0]; % yellow color state for death or unpeopled
vc(2,:) = [0,1,1]; % black color state for curative
vc(3,:) = [0,0,1]; % blue color state for delitescent
vc(4,:) = [1,0,0]; % red color state for catch an illness
vc(5,:) = [0,1,0]; % green color state for normal people
vc(6,:) = [1,0,1]; % black color state for normal people
% Sc = ones(2N+1,2N+1,3); % color state matrix
Sc = mark_color(S,vc);
%%%%%% initialization of Figure
set(gcf,‘Position’,[200,60,600,700]);
axes(‘Position’,[0.1,0.4,0.6,0.5]); % mapping
I1 = imshow(Sc);
th = title(‘The spread days of the disease: 1’,‘FontSize’,14);
axes(‘Position’,[0.78,0.4,0.2,0.5]); hold on;
ylim([-0.1,1.05]); xlim([-0.1,1]);
% box on;
axis off
rectangle(‘Position’,[0,0.35,0.3,0.10],‘FaceColor’,vc(1,:));
text(0.4,0.4,{‘D’});
rectangle(‘Position’,[0,0.05,0.3,0.10],‘FaceColor’,vc(2,:));
text(0.4,0.1,{‘R’});
rectangle(‘Position’,[0,0.65,0.3,0.10],‘FaceColor’,vc(3,:));
text(0.4,0.70,{‘E’});
rectangle(‘Position’,[0,0.5,0.3,0.10],‘FaceColor’,vc(4,:));
text(0.4,0.55,{‘I’});
rectangle(‘Position’,[0,0.8,0.3,0.10],‘FaceColor’,vc(5,:));
text(0.4,0.85,{‘S’});
rectangle(‘Position’,[0,0.2,0.3,0.10],‘FaceColor’,vc(6,:));
text(0.4,0.25,{‘Q’});
set(gca,‘Xtick’,[],‘Ytick’,[]);
axes(‘Position’,[0.1,0.08,0.6,0.2]); % stats
Bh(1) = bar(1,sum(sum(S1)),‘FaceColor’,vc(2,:)); hold on;
Bh(2) = bar(2,sum(sum(S2)),‘FaceColor’,vc(3,:));
Bh(3) = bar(3,sum(sum(S==3)),‘FaceColor’,vc(4,:));
ylabel(‘number of people’);
xlim([0,4]);set(gca,‘Xtick’,[0:4])
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% iterative calculation
S(N+N/2,N+1) = 3; % center cell is catched an illness
S(N+(N-10)/2,N+1) = 3; % center cell is catched an illness
S(N-N/2,N+1) = 3; % center cell is catched an illness
S(N+1,N+1) = 3; % center cell is catched an illness
S(N+1,N+N/2) = 3; % center cell is catched an illness
S(N+1,N-N/2) = 3; % center cell is catched an illness
SL(S3) = 0; % time for the state of catching an illness
SGL = zeros(2*N+1); % 隔離時間
Sc = mark_color(S,vc);
set(I1,‘CData’,Sc);
Yao = 60;
cure=0;
ge=0;
ill=0;
for ktime = 1:200; % unit: day
% 潛伏者 5-10天有 80%的概率變成發(fā)病者,也就是狀態(tài)2有80%概率變成狀態(tài)3 ,
% 有20%的概率變成狀態(tài)1
s2 = find(S2);
for k = 1:length(s2);
if SL(s2(k))>4.5 & SL(s2(k))<10.5; % 5-10天
if rand<0.8;
S(s2(k)) = 3; % 2–>3, 80%
ill=ill+1;
SL(s2(k)) = 1; % time updata
else
S(s2(k)) = 1; % 2–>1, 20%
cure=cure+1;%治好人數(shù)更新
SL(s2(k)) = 1; % time updata
end
else
SL(s2(k)) = SL(s2(k))+1; % updating data
end
end
s3 = find(S3);
%%隔離死亡
s5 = find(S5);
for k = 1:length(s5);
if SGL(s5(k))>0.5 & SGL(s5(k))<10.5; % 5-10天
if rand<0.4;
S(s5(k)) = 0; % 2–>0, 80%
Dn = Dn +1;
% ill=ill+1;
SGL(s5(k)) = 0; % time updata
else
% S(s5(k)) = 1; % 2–>1, 20%
SGL(s5(k)) = SGL(s5(k))+1; % time updata
end
else
if rand<0.5;
S(s5(k)) = 0; % 2–>3, 80%
Dn = Dn +1;
% ill=ill+1;
SGL(s5(k)) = 0; % time updata
else
% S(s5(k)) = 1; % 2–>1, 20%
SGL(s5(k)) = SGL(s5(k))+1; % time updata
end
% SL(s2(k)) = SL(s2(k))+1; % updating data
end
end
%%%隔離狀態(tài)
for k = 1:length(s3);
if ktime>60;
if rand<min(SL(s3(k))*0.05+0.05,0.4);
S(s3(k)) = 5;
ge = ge+1;
if SL(s3(k))<10.5 % 如果在患病 10天內(nèi)發(fā)現(xiàn)
if rand<0.8 & rand<min(1,Yao/ge)
S(s3(k)) = 1; % 5–>1. % 50%被治好, 藥物的期望0.9,也就是有90%可以被治好
cure=cure+1; %治好人數(shù)更新
SGL(s3(k)) = 0; % 隔離時間清空
else
SGL(s3(k)) = SGL(s3(k))+1; % 隔離時間更新
end
elseif SL(s3(k))>10.5 % 患病10天到20天
if rand<0.2 & rand<rand<min(1,Yao/ge);
S(s3(k)) = 1; % 5–>1. 有20%被治好,藥物的期望0.9,也就是有90%可以被治好
cure=cure+1;%治好人數(shù)更新
SGL(s3(k)) = 0; % 隔離時間清空
else
SGL(s3(k)) = SGL(s3(k))+1; % 隔離時間更新
end
elseif SGL(s3(k))<10.5 % 在隔離5天內(nèi)
if rand<0.2 & rand<rand<min(1,Yao/ge)
S(s3(k)) = 1; % 5–>1. 有20%被治好,藥物的期望0.9,也就是有90%可以被治好
cure=cure+1;%治好人數(shù)更新
SGL(s3(k)) = 0; % 隔離時間清空
else
SGL(s3(k)) = SGL(s3(k))+1; % 隔離時間更新
end
elseif SGL(s3(k))>10.5 % 在隔離5天以上
if rand<0.5 & rand<rand<min(1,Yao/ge)
S(s3(k)) = 1; % 5–>1. 有50%被治好,藥物的期望0.9,也就是有90%可以被治好
cure=cure+1;%治好人數(shù)更新
SGL(s3(k)) = 0; % 隔離時間清空
else
% S(s3(k)) = 0; % 5–>1. 有50%未被治好,藥物沒有作用而死亡
% SGL(s3(k)) = 0; % 隔離時間清空
SGL(s3(k)) = SGL(s3(k))+1; % 隔離時間更新
end
end
else
SL(s3(k)) = SL(s3(k))+1; % 未隔離者發(fā)表時間更新
end
end
三、運(yùn)行結(jié)果
四、matlab版本及參考文獻(xiàn)
1 matlab版本
2014a
2 參考文獻(xiàn)
[1] 包子陽,余繼周,楊杉.智能優(yōu)化算法及其MATLAB實(shí)例(第2版)[M].電子工業(yè)出版社,2016.
[2]張巖,吳水根.MATLAB優(yōu)化算法源代碼[M].清華大學(xué)出版社,2017.
[3]【數(shù)學(xué)建模】元胞自動機(jī).博主:二進(jìn)制 人工智能
總結(jié)
以上是生活随笔為你收集整理的【元胞自动机】基于matlab元胞自动机模拟SEIR传播模型【含Matlab源码 2156期】的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 处理机的三级调度
- 下一篇: 计算机二级access选择题技巧,计算机