logistic模型及其matlab算法

1、Logistic曲线的三种参数估计方法作者QQ：2377389590Logistic曲线的参数估计1844或1845年，比利时数学家PierreFrangoisVerhulst提出了logistic方程，这是一个对S型曲线进行数学描述的模型。一百多年来，这个方程多次应用于一些特殊的领域建模与预测，例如单位面积内某种生物的数量、人口数量等社会经济指标、某种商品(例如手机)的普及率等。logistic方程定义如下：1x(1)tc+aebt其中，t通常表示时间变量，a,b和c为模型的参数；当趋势比较完整时a0b0,c,0o其曲线如图1所示：根据1图和(1)方程式得：当t-时，x(t)0；当t+时，x

2、1/co为了研究Logistic曲线的增长特性，对(1)式求导一阶导数得：dx-abebt小,0dt(c+aebt)2图1logistic方程的曲线示意图设xt(t=1,2,，n)为观测样本，对于logistic方程的参数，常规的估计方法有三种。Yule算法：根据式(1)，有(7)t1xt1xt1caeb(t1)(2)caebt(aebt+c-c)(1-eb)caebt，(1-eb)-c(1-eb)xtx-xc设z二占F、丫，1-eb以及，-c（1-eb），则式变形为线性方程z，Y+x,txttt1利用普通最小二乘（OLS）方法可以得到这个方程参数的估计值，b和c的估计值也可以进一步得到。为得

3、到a的估计值，将式（1）变形为：ln-c，ina+bt(二丿(3)(7)(7)左右分别对t求和1)-c，nIna+t=11xt丿a二exp厂1in1)-cn一t,11xt丿因此，a的估计值为：n(n1)b2n(n+1)b2(4)(5)(7)(7)Rhodes算法根据式（1），有1，c+aeb(t1)xt1，c-ceb+ceb+aeb(t+1)(6)eb，c(1-eb)+一x1设z=一txt11s，、丫，txttc（1-eb）以及，eb，则式（6）变形为z，丫+s。用普tt通最小二乘（OLS）方法得到这个方程参数的估计值，并进一步得到b和c的估计值，然后利用式（3）-式（5）的方法得到a的估计值

4、。Nair算法：式（2）可以进一步写为：i+4_二(1-eb)(1-cx)xtt1(8)即：1x(1-cx)=t1t HYPERLINK l bookmark361一ebx一xt1t因此，有1+eb22x(1-cx)，+4t1xxt1tx2cxx+xt+1t+1tt(9)t+1t+1t+1式（9）可以整理为：111eb2c(1eb)1+eb(10)t+11111设z=、s，+txxtxxtt+1tt+1xxtt+12c(1eb)1ebr，以及卩，则式(10)可写1+eb1+eb(8)(8)为z，丫+卩s。用普通最小二乘（OLS）方法得到这个方程参数的估计值，并进一步得到ttb和c的估计值，然后

5、利用式（3）-式（5）的方法得到a的估计值。2.列题：中国1965-2011年CO2排放量（单位：亿吨）如表1所示:表1中国1965-2011年CO2排放量1965196619671968196919701971197219731974480.9522.0468.8469.5573.8737.8869.8933.7977.2997.719751976197719781979198019811982198319841120.31176.11284.81422.11462.11499.71473.11539.21637.01771.019851986198719881989199019911992

6、199319941886.51994.62145.72292.02396.82387.02484.42580.82750.22915.719951996199719981999200020012002200320043163.83231.93319.53319.63484.03550.63613.93833.14471.25283.020052006200720082009201020115803.26415.56797.97033.57636.38209.88979.1用logistic方程模拟我国CO2排放量的变化趋势，分别用三种方法估计方程参数，并分别计算三种方法的MAPE及未来五年CO

7、2排放量的预测结果。2.1.Yule算法：clear;clc;%Yule算法：%author:朱伟杰%date:2018-1-24X=480.9,522,468.8,469.5,573.8,737.8,869.8,933.7,977.2,.3550.6,3613.9,3833.1,4471.2,5283,.5803.2,6415.5,6797.9,7033.5,7636.3,8209.8,8979.1n=length(X)-1fort=1:nZ(t)=(X(t+1)-X(t)/X(t+1)endX1=ones(46,1)X(1:n)Y=ZB,Bint,r,rint,stats=regress(

8、Y,Xl)%最小二乘(0LS)gamma=B(1,1)beta=B(2,1)b=log(1-gamma)c=beta/(exp(b)-1)a=exp(sum(log(1./X(1:n)-c)-n*(n+1)*b/2)/n)XX=1965:2016YY=1./(c+a*exp(b*(XX-1965)plot(XX,YY,r-o)holdonplot(XX(1:length(X),X,g-入)legend(预测值,实际值)xlabel(年份);ylabel(C0_2排放量);title(C0_2预测值和实际值曲线图(Yule法)set(gca,XTick,1965:2:2017)gridonfor

9、matshort;forecast=YY(end-4:end)%C02排放量的预测结果MAPE=sum(abs(YY(1:n+1)-X)./X)/length(X)%平均相对差值a,b,c2.2.Rhodes算法：clear;clc;%Rhodes算法%author:朱伟杰%date:2018-1-24X=480.9,522,468.8,469.5,573.8,737.8,869.8,933.7,977.2,.3550.6,3613.9,3833.1,4471.2,5283,.n=length(X)-1fort=1:nZ(t)=1/X(t+1)S(t)=1/X(t)endX1=ones(46,

10、1)S(1:n)Y=ZB,Bint,r,rint,stats=regress(Y,Xl)%最小二乘(OLS)gamma=B(1,1)beta=B(2,1)b=log(beta)c=gamma/(1-exp(b)a=exp(sum(log(1./X(1:n+1)-c)-(n+1)*(n+2)*b/2)/(n+1)XX=1965:2016YY=1./(c+a*exp(b*(XX-1965)plot(XX,YY,r-o)holdonplot(XX(1:length(X),X,k一入)set(gca,XTick,1965:2:2017)legend(预测值,实际值)xlabel(年份);ylabel(

11、CO_2排放量);title(C0_2预测值和实际值曲线图(Rhodes法)gridonformatshort;forecast=YY(end-4:end)%C02排放量的预测结果MAPE=sum(abs(YY(1:n+1)-X)./X)/length(X)%平均相对差值a,b,c匚丐预测也和实麻由曲线图很hod歸袪jmT1辨13711旳519771W11W1*9E忖怡11卿?S的阳2WJ?W5旳OS20n却T左栏.8Nair算法：clear;clc;%Nair算法%author:朱伟杰%date:2018-1-24X=480.9,522,468.8,469.5,573.8,737.8,869

12、.8,933.7,977.2,.TOC o 1-5 h z.3550.6.3613.9.3833.1.4471.2.5283.n=length(X)-1fort=1:nZ(t)=1/X(t)-1/X(t+1)S(t)=1/X(t)+1/X(t+1)endX1=ones(46,1)S(1:n)Y=ZB,Bint,r,rint,stats=regress(Y,X1)%最小二乘(OLS)gamma=B(1,1)beta=B(2,1)b=log(1-beta)/(1+beta)c=gamma*(1+exp(b)/(2*(exp(b)-1)a=exp(sum(log(1./X(1:n)-c)-n*(n+

13、1)*b/2)/n)XX=1965:2016YY=1./(c+a*exp(b*(XX-1965)plot(XX,YY,r-o)holdonplot(XX(1:length(X),X,g-入)legend(预测值,实际值)xlabel(年份);ylabel(二氧化碳排放量);title(二氧化碳预测值和实际值曲线图(Nair法)set(gca,XTick,1965:2:2017)gridonformatshort;forecast=YY(end-4:end)%C02排放量的预测结果MAPE=sum(abs(YY(1:n+1)-X)./X)/length(X)%平均相对差值a,b,c1WSS1!M|IMSIWrlIMiI-IJ.S-JJirIJUUb?Ett-d5WIWiWJ1J/J111TJ21旧匕IJ-J7TJLSh：I髭11ZD!-2Ph*一軌化轨预敌恒和实际恒曲绘图小謝惟）1BKD5ZW-预浏悄实际勺I根据编程计算，三种方法参数估计的结果及MAPE如表1所示：表1三种方法的参数估计结果及误差分析结果abcMAPEYule算法0.0020-0.0580-0.0000245790.1141Rhodes算法0.0019-0.0802