产品
听说大家都挺无聊的，我也是，于是昨天晚上我打算自己动手用matlab算一下我啥时候能结束这段无聊的宅家生活。
宅家之前，虽然也“无聊”但是我可以骑着我心爱的小电动在这个小城里左右兜转，有时候还会去拍拍照片，喝喝奶茶，和朋友聚聚餐，看看电影啥的o(*￣▽￣*)ブ，但是现在在家就是卧室，客厅，厕所三点一线的生活。/(ㄒoㄒ)/~~
宅家之后，只能通过王者峡谷和各位进行沟通互动，虽然也玩过弹弹堂，洛克王国，和平精英，但是主要方式还是王者农药。
樟树市启动一级响应的第一天，封闭生活算是正式开始了，本来人满为患的时代广场（樟树市的时代广场）转眼变成空空如也，可以说全国都进入了紧张阶段，武汉人民的情况更是牵动着全国人们的心。
一级响应前几天：
时代广场
时代广场
一级响应启动后：
看到冷清的街头我常常想这是怎么了，为什么短短几天世界变成了我不认识的样子，一切就好像在放电影一样，和那位抽泣着说，”我的武汉没有车了“ 的记者一样，看到大过年冷清的樟树街头，我就不由得心里就泛酸。
最近几天开始想找事情干，录了几个游戏视频，学pr剪视频，用matlab图像处理玩了一下当年红极一时的bad apple，反正都是以前没玩过的东西。
昨天看到一个关于matlab预测疫情的视频，于是想自己也动手做一下，虽然网上已经有很多大佬分析过了各种数据网上也能查得到但是还是忍不住自己动手，第一是实在想找事情做而且是要自己感兴趣的，第二是想学习一下数学建模的过程（本人菜鸡单纯娱乐一下）。
2月17日晚上，我查找了国家卫健委发布的从1月20日到2月16日疫情通报的所有数据，并且进行了部分统计作为准备工作。
来源：国家卫健委
参考了b站这位up的建模视频
以及毕导的这个视频
最后查阅了知网的相关论文（感谢江苏大学与知网合作让我能免费下载论文），主要是下面这一篇。
最后用matlab动手操作，通过函数拟合，我惊讶的发现康复率成幂函数增长，最近更超过了1.0%，死亡率成线性增长，这不就说明疫情会变好的越来越快吗。
o(*￣▽￣*)ブ
康复率
死亡率
根据最近的一些情况对SEIR模型进行同步修正，比如2月2日国家进行管控后，感染者和潜伏着的接触人数大大减少，现在的治疗效率大大提高，康复率利用幂函数同步调整而不是原来的单调的线性增长，当然建模的一些数据还是比较粗糙的，但是实际最后结果只会比我的模型结果更好不会差。代码如下感兴趣的朋友可以自己去试一试。
这是建模代码：
clear;clc
S=8837300/8; %易感人群
E=2.4; %潜伏期人群
I=1; %感染者
R=0; %康复者
N=S+E+I+R; %总数
r1=20; %潜伏期患者每天接触到的人
r2=10; %感染者每天接触到的人
b1=0.045; %潜伏期患者传染概率
b2=0.045; %感染者传染概率
a=0.10; %潜伏期患者变成感染者
%r患者康复概率
%d患者死亡概率
[rr,d] = getData_COVID_19( );
days=200;
for i=1:200
if i>=16
r1=3;
r2=5;
end
if i<=28
r=rr(i);
else
r=rr(28);
end
S(i+1)=S(i)-(r1*b1*E(i)*(S(i)/N)+r2*b2*I(i)*(S(i)/N));
E(i+1)=E(i)+(r1*b1*E(i)*(S(i)/N)+r2*b2*I(i)*(S(i)/N))-(E(i)*a);
I(i+1)=I(i)+E(i)*a-(I(i)*r+I(i)*d);
R(i+1)=R(i)+I(i)*r;
N=S(i+1)+E(i+1)+ I(i+1)+R(i+1);
end
x=1:days+1;
subplot(2,1,1)
plot(x,S,x,E,x,I,x,R);
xlabel('天数');
ylabel('人数');
legend('易感人群','潜伏期患者','感染者','康复者');
hold on;
subplot(2,1,2)
plot(rr(1:28))
xlabel('天数');
ylabel('康复率');
这是康复率函数拟合代码：
function [fitresult, gof] = get_RE_2(pat, pat_re)
%CREATEFIT(PAT,PAT_RE)
% Create a fit.
%
% Data for 'untitled fit 1' fit:
% X Input : pat
% Y Output: pat_re
% Output:
% fitresult : a fit object representing the fit.
% gof : structure with goodness-of fit info.
% 另请参阅 FIT, CFIT, SFIT.
% 由 MATLAB 于 18-Feb-2020 11:28:16 自动生成
%% Fit: 'untitled fit 1'.
[xData, yData] = prepareCurveData( pat, pat_re );
% Set up fittype and options.
ft = fittype( 'power1' );
opts = fitoptions( 'Method', 'NonlinearLeastSquares' );
opts.Display = 'Off';
opts.StartPoint = [0.162869601184208 0.900104485756097];
% Fit model to data.
[fitresult, gof] = fit( xData, yData, ft, opts );
这是死亡率线性拟合代码：
function [fitresult, gof] = getDIE(pat, pat_die)
%CREATEFIT(PAT,PAT_DIE)
% Create a fit.
%
% Data for 'untitled fit 1' fit:
% X Input : pat
% Y Output: pat_die
% Output:
% fitresult : a fit object representing the fit.
% gof : structure with goodness-of fit info.
% 另请参阅 FIT, CFIT, SFIT.
% 由 MATLAB 于 17-Feb-2020 20:33:33 自动生成
%% Fit: 'untitled fit 1'.
[xData, yData] = prepareCurveData( pat, pat_die );
% Set up fittype and options.
ft = fittype( 'poly1' );
% Fit model to data.
[fitresult, gof] = fit( xData, yData, ft );
这是同步更新死亡率和康复率代码：
function [rr,d] = getData_COVID_19( )
[data,~,~]=xlsread('.\data\全国新冠肺炎数据.xls');
pat=data(:,2); %确诊人数
pat_re=data(:,5); %康复人数
pat_die=data(:,4); %死亡人数
[fitresult, gof] = get_RE_2(pat, pat_re);
m=fitresult.a;
n=fitresult.b;
for i=1:length(pat)
rr(i)=m*pat(i)^(n-1);
end
[fitresult, gof] = getDIE(pat, pat_die);
d=fitresult.p1;
end
最后matlab运行结果如下图，当然如果想自己搞丁香医生里面那些折线图只要需要在我的代码和数据里直接plot画出来就好了，对我个人没啥意义我就没画了。
疫情走势图
通过分析上图可以发现54天也就是3月10日左右感染人数将达到峰值，最近国家调控越来越严格，火神山，雷神山等医院建成，检验试剂的研发以及ct检查的应用，实际到来日期肯定会在3月10日之前，最终人数也会远远小于我模拟的14万，个人保守估计是10万人以内。现在我只希望越来越好，等到春暖花开那一天，一切回到什么也没发生过那时候的状态，你上你的班，我上我的课。
最后关键时刻，大家一定要耐得住寂寞，守住我们最后的胜利，病毒散去！武汉必胜！人类必胜！最后，致敬日日夜夜工作在一线的医护人员，你们是我心中的英雄，希望每位白衣天使都能被温柔以待。
武汉加油
155-2731-8020
