• 数学建模--预测类模型


    目录

    一、中短期预测

    1、灰色预测法

    ①适用范围

    ②模型实现

     2、回归分析

    ①适用范围

    ②模型实现

     3、时间序列分析

    ①自适应滤波法

    ②指数平滑法

    ③移动平均法

    4、微分方程

    二、长期预测

    1、神经网络预测

    2、logistic模型

    ①模型介绍

    ②模型分析及代码

    一、中短期预测

    1、灰色预测

    灰色预测模型 ( Gray Forecast Model )是通过少量的、不完全的信息,建立数学模型并做出预测的一种预测方法。当我们应用运筹学的思想方法解决实际问题,制定发展战略和政策、进行重大问题的决策时,都必须对未来进行科学的预测 。预测是根据客观事物的过去和现在的发展规律,借助于科学的方法对其未来的发展趋势和状况进行描述和分析,并形成科学的假设和判断。

    ①适用范围

    适用范围:该模型使用的不是原始数据的序列,而是生成的数据序列。核心体系是Grey Model,即对原始数据作累加生成(或其他处理生成)得到近似的指数规律再进行建模的方法。

    优点:在处理较少的特征值数据,不需要数据的样本空间足够大,就能解决历史数据少、序列的完整性以及可靠性低的问题,能将无规律的原始数据进行生成得到规律较强的生成序列。

    缺点:只适用于中短期的预测,只适合近似于指数增长的预测。

    ②模型实现

    GM(1,1)模型适用 于具有较强指数规律的序列,只能描述单调的变化过程,而 Verhulst 模型则适用于非单 调的摆动发展序列或具有饱和状态的 S 形序列。下面对比一下两种方法:

    用matlab语言采用Verhulst模型预测

    1. clc,clear
    2. x1=[4.93 5.33 5.87 6.35 6.63 7.15 7.37 7.39 7.81 8.35 9.39 10.59 10.94 10.44];
    3. n=length(x1);
    4. nian=1990:2003;
    5. plot(nian,x1,'o-');
    6. x0=diff(x1);
    7. x0=[x1(1),x0]
    8. for i=2:n
    9. z1(i)=0.5*(x1(i)+x1(i-1));
    10. end
    11. z1
    12. B=[-z1(2:end)',z1(2:end)'.^2]
    13. Y=x0(2:end)'
    14. abhat=B\Y %估计参数 a,b 的值
    15. x=dsolve('Dx+a*x=b*x^2','x(0)=x0'); %求解常微分方程
    16. x=subs(x,{'a','b','x0'},{abhat(1),abhat(2),x1(1)}); %代入参数值
    17. yuce=subs(x,'t',0:14) %计算预测值
    18. digits(6); x=vpa(x) %显示微分方程的解,为了提高计算精度,把该语句放在计算预测值之后,或者不使用该语句
    19. yuce(16)=yuce(15);
    20. x1_all=[x1,9.92,10.71];
    21. epsilon=x1_all-yuce %计算残差
    22. delta=abs(epsilon./x1_all) %计算相对误差
    23. delta_mean=mean(delta) %计算平均相对误差
    24. x1_all_0=x1_all-x1_all(1); %数据列的始点零化像
    25. yuce_0=yuce-yuce(1); %数据列的始点零化像
    26. s0=abs(sum(x1_all_0(1:end-1))+0.5*x1_all_0(end));
    27. s1=abs(sum(yuce_0(1:end-1))+0.5*yuce_0(end));
    28. tt=yuce_0-x1_all_0;
    29. s1_s0=abs(sum(tt(1:end-1))+0.5*tt(end));
    30. absdegree=(1+s0+s1)/(1+s0+s1+s1_s0) %计算灰色绝对关联度
    31. c=std(epsilon,1)/std(x1_all,1) %计算标准差比值

    参考文章:http://t.csdn.cn/OIMNK参考文章:http://t.csdn.cn/OIMNK

    用python语言采用GM(1,1)微分方程预测

    1. import numpy as np
    2. import math
    3. history_data = [4.93, 5.33, 5.87, 6.35, 6.63, 7.15, 7.37, 7.39, 7.81, 8.35,
    4. 9.39, 10.59, 10.94, 10.44]
    5. n = len(history_data)
    6. X0 = np.array(history_data)
    7. # 累加生成
    8. history_data_agg = [sum(history_data[0:i + 1]) for i in range(n)]
    9. X1 = np.array(history_data_agg)
    10. # 计算数据矩阵B和数据向量Y
    11. B = np.zeros([n - 1, 2])
    12. Y = np.zeros([n - 1, 1])
    13. for i in range(0, n - 1):
    14. B[i][0] = -0.5 * (X1[i] + X1[i + 1])
    15. B[i][1] = 1
    16. Y[i][0] = X0[i + 1]
    17. # 计算GM(1,1)微分方程的参数a和u
    18. # A = np.zeros([2,1])
    19. A = np.linalg.inv(B.T.dot(B)).dot(B.T).dot(Y)
    20. a = A[0][0]
    21. u = A[1][0]
    22. # 建立灰色预测模型
    23. XX0 = np.zeros(n)
    24. XX0[0] = X0[0]
    25. for i in range(1, n):
    26. XX0[i] = (X0[0] - u / a) * (1 - math.exp(a)) * math.exp(-a * (i));
    27. # 模型精度的后验差检验
    28. e = 0 # 求残差平均值
    29. for i in range(0, n):
    30. e += (X0[i] - XX0[i])
    31. e /= n
    32. # 求历史数据平均值
    33. aver = 0;
    34. for i in range(0, n):
    35. aver += X0[i]
    36. aver /= n
    37. # 求历史数据方差
    38. s12 = 0;
    39. for i in range(0, n):
    40. s12 += (X0[i] - aver) ** 2;
    41. s12 /= n
    42. # 求残差方差
    43. s22 = 0;
    44. for i in range(0, n):
    45. s22 += ((X0[i] - XX0[i]) - e) ** 2;
    46. s22 /= n
    47. # 求后验差比值
    48. C = s22 / s12
    49. # 求小误差概率
    50. cout = 0
    51. for i in range(0, n):
    52. if abs((X0[i] - XX0[i]) - e) < 0.6754 * math.sqrt(s12):
    53. cout = cout + 1
    54. else:
    55. cout = cout
    56. P = cout / n
    57. if (C < 0.35 and P > 0.95):
    58. # 预测精度为一级
    59. m = 2 # 请输入需要预测的年数
    60. # print('往后m各年负荷为:')
    61. f = np.zeros(m)
    62. for i in range(0, m):
    63. f[i] = (X0[0] - u / a) * (1 - math.exp(a)) * math.exp(-a * (i + n))
    64. print(f[i])
    65. else:
    66. print('灰色预测法不适用')

     2、回归分析

    ①适用范围

         回归分析是研究变量之间因果关系的一种统计模型;因变量就是结果,自变量就是原因;基于结果变量(因变量)的种类,回归分析可分为:线性回归(因变量为连续变量)、logistic回归(因变量为分类变量)、柏松回归(因变量为计数变量);这三种回归模型中自变量则可以是任意类型的变量。

      有的自变量对因变量的影响不是很大,且自变量之间可能存在多重共线性(即可能不完全独立),通过建立逐步回归分析,进行X因子筛选。

    比如:收入水平于受教育程度、所在行业、工作年限、工作种类的关系;公路客运量与人口增长量、私家车保有量、国民生产总值、国民收入等因素的关系。

    ②模型实现

    几种回归分析的matlab代码如下

    非线性回归

    1. %% 非线性回归
    2. clc, clear
    3. x0 = [1 8.55 470 300 10
    4. 2 3.79 285 80 10
    5. 3 4.82 470 300 120
    6. 4 0.02 470 300 120
    7. 5 2.75 470 80 10
    8. 6 14.39 100 190 10
    9. 7 2.54 100 80 65
    10. 8 4.35 470 190 65
    11. 9 13 100 300 54
    12. 10 8.5 100 300 120
    13. 11 0.05 100 80 120
    14. 12 11.32 285 300 10
    15. 13 3.13 285 190 120];
    16. x = x0(:, 3:5);
    17. y = x0(:, 2);
    18. % 参数的初始预估值
    19. beta = [0.1, 0.05, 0.02, 1, 2]';
    20. [betahat, r, j] = nlinfit(x, y, @func, beta);
    21. betaci = nlparci(betahat, r, 'jacobian', j);
    22. % 回归系数以及置信区间
    23. betaa = [betahat, betaci];
    24. % y的预测值以及置信区间半径
    25. [yhat, delta] = nlpredci(@func, x, betahat, r, 'jacobian', j);
    26. % 绘制交互式画面
    27. nlintool(x, y, 'func', beta)

     逐步回归模型

    1. %% 逐步回归模型
    2. clc, clear
    3. x0 = [1 8.55 470 300 10
    4. 2 3.79 285 80 10
    5. 3 4.82 470 300 120
    6. 4 0.02 470 300 120
    7. 5 2.75 470 80 10
    8. 6 14.39 100 190 10
    9. 7 2.54 100 80 65
    10. 8 4.35 470 190 65
    11. 9 13 100 300 54
    12. 10 8.5 100 300 120
    13. 11 0.05 100 80 120
    14. 12 11.32 285 300 10
    15. 13 3.13 285 190 120];
    16. x = x0(:, 3:5);
    17. y = x0(:, 2);
    18. stepwise(x, y)

     一元多项回归模型

    1. %% 一元多项式回归模型
    2. % 选取 y = a2*x^2 + a1*x + a0
    3. clc, clear
    4. x0 = 17: 2: 29;
    5. y0 = [20.48, 25.13, 26.15, 30.0, 26.1, 20.3, 19.35];
    6. [p, s] = polyfit(x0, y0, 2); % 2表示阶数
    7. p % p返回a2,a1,a0三个参数
    8. % 求polyfit所得的回归多项式在x0处的预测值Y
    9. %及预测值的显著性为1-alpha的置信区间DELTA;
    10. %alpha缺省时为0.05它假设polyfit函数数据输入的误差是独立正态的,并且方差为常数
    11. [y, delta] = polyconf(p, x0, s); % delta为置信区间半径
    12. polytool(x0, y0, 2)

     

     多元二项式回归

    1. %% 多元二项式回归
    2. % rstool(x, y, model, alpha)
    3. % model有四种:线性,纯二次,交叉,完全二次
    4. clc, clear
    5. x1 = [120 140 190 130 155 175 125 145 180 150]';
    6. x2 = [100 110 90 150 210 150 250 270 300 250]';
    7. y = [102 100 120 77 46 93 26 69 65 85]';
    8. x = [x1 x2];
    9. rstool(x, y, 'purequadratic')

    一元(多元)线性回归

    1. %% 一元(多元)线性回归
    2. %x = input('输入需要预测的一组数据:(如[0.10 0.11 0.12 ...)\n');
    3. % 测试用例:
    4. x1 = [0.1:0.01:0.18]';
    5. %y = input('输入需要预测的一组数据:(如[0.10 0.11 0.12 ...)\n');
    6. % 测试用例:
    7. y = [42, 41.5, 45.0, 45.5, 45.0, 47.5, 49.0, 55.0, 50.0]';
    8. plot(x1, y, '+') % 观察是否近似为线性模型
    9. x = [ones(9, 1), x1];
    10. [b, bint, r, rint, stats] = regress(y, x);
    11. b %回归系数估计值
    12. bint %回归系数置信区间
    13. stats %检验回归系统的统计量:复判定系数等
    14. rcoplot(r, rint) %残差向量与残差置信区间
    15. %--根据rcoplot的绘制结果,剔除异常点之后重新生成回归系统--%

     

    偏最小二乘回归(PLS)   python代码

    1. import numpy as np
    2. x1 = [191, 189, 193, 162, 189, 182, 211, 167, 176, 154, 169, 166, 154, 247, 193, 202, 176, 157, 156, 138]
    3. x2 = [36, 37, 38, 35, 35, 36, 38, 34, 31, 33, 34, 33, 34, 46, 36, 37, 37, 32, 33, 33]
    4. x3 = [50, 52, 58, 62, 46, 56, 56, 60, 74, 56, 50, 52, 64, 50, 46, 62, 54, 52, 54, 68]
    5. y1 = [5, 2, 12, 12, 13, 4, 8, 6, 15, 17, 17, 13, 14, 1, 6, 12, 4, 11, 15, 2]
    6. y2 = [162, 110, 101, 105, 155, 101, 101, 125, 200, 251, 120, 210, 215, 50, 70, 210, 60, 230, 225, 110]
    7. y3 = [60, 60, 101, 37, 58, 42, 38, 40, 40, 250, 38, 115, 105, 50, 31, 120, 25, 80, 73, 43]
    8. # -----数据读取
    9. data_raw = np.array([x1, x2, x3, y1, y2, y3])
    10. data_raw = data_raw.T # 输入原始数据,行数为样本数,列数为特征数
    11. # -----数据标准化
    12. num = np.size(data_raw, 0) # 样本个数
    13. mu = np.mean(data_raw, axis=0) # 按列求均值
    14. sig = (np.std(data_raw, axis=0)) # 按列求标准差
    15. data = (data_raw - mu) / sig # 标准化,按列减去均值除以标准差
    16. # -----提取自变量和因变量数据
    17. n = 3 # 自变量个数
    18. m = 3 # 因变量个数
    19. x0 = data_raw[:, 0:n] # 原始的自变量数据
    20. y0 = data_raw[:, n:n + m] # 原始的变量数据
    21. e0 = data[:, 0:n] # 标准化后的自变量数据
    22. f0 = data[:, n:n + m] # 标准化后的因变量数据
    23. # -----相关矩阵初始化
    24. chg = np.eye(n) # w到w*变换矩阵的初始化
    25. w = np.empty((n, 0)) # 初始化投影轴矩阵
    26. w_star = np.empty((n, 0)) # w*矩阵初始化
    27. t = np.empty((num, 0)) # 得分矩阵初始化
    28. ss = np.empty(0) # 或者ss=[],误差平方和初始化
    29. press = [] # 预测误差平方和初始化
    30. Q_h2 = np.zeros(n) # 有效性判断条件值初始化
    31. # -----求解主成分
    32. for i in range(n): # 主成分的总个数小于等于自变量个数
    33. # -----求解自变量的最大投影w和第一主成分t
    34. matrix = e0.T @ f0 @ f0.T @ e0 # 构造矩阵E'FF'E
    35. val, vec = np.linalg.eig(matrix) # 计算特征值和特征向量
    36. index = np.argsort(val)[::-1] # 获取特征值从大到小排序前的索引
    37. val_sort = val[index] # 特征值由大到小排序
    38. vec_sort = vec[:, index] # 特征向量按照特征值的顺序排列
    39. w = np.append(w, vec_sort[:, 0][:, np.newaxis], axis=1) # 储存最大特征向量
    40. w_star = np.append(w_star, chg @ w[:, i][:, np.newaxis], axis=1) # 计算 w*的取值
    41. t = np.append(t, e0 @ w[:, i][:, np.newaxis], axis=1) # 计算投影
    42. alpha = e0.T @ t[:, i][:, np.newaxis] / (t[:, i] @ t[:, i]) # 计算自变量和主成分之间的回归系数
    43. chg = chg @ (np.eye(n) - (w[:, i][:, np.newaxis] @ alpha.T)) # 计算 w 到 w*的变换矩阵
    44. e1 = e0 - t[:, i][:, np.newaxis] @ alpha.T # 计算残差矩阵
    45. e0 = e1 # 更新残差矩阵
    46. # -----求解误差平方和ss
    47. beta = np.linalg.pinv(t) @ f0 # 求回归方程的系数,数据标准化,没有常数项
    48. res = np.array(f0 - t @ beta) # 求残差
    49. ss = np.append(ss, np.sum(res ** 2)) # 残差平方和
    50. # -----求解残差平方和press
    51. press_i = [] # 初始化误差平方和矩阵
    52. for j in range(num):
    53. t_inter = t[:, 0:i + 1]
    54. f_inter = f0
    55. t_inter_del = t_inter[j, :] # 把舍去的第 j 个样本点保存起来,自变量
    56. f_inter_del = f_inter[j, :] # 把舍去的第 j 个样本点保存起来,因变量
    57. t_inter = np.delete(t_inter, j, axis=0) # 删除自变量第 j 个观测值
    58. f_inter = np.delete(f_inter, j, axis=0) # 删除因变量第 j 个观测值
    59. t_inter = np.append(t_inter, np.ones((num - 1, 1)), axis=1)
    60. beta1 = np.linalg.pinv(t_inter) @ f_inter # 求回归分析的系数,这里带有常数项
    61. res = f_inter_del - t_inter_del[:, np.newaxis].T @ beta1[0:len(beta1) - 1, :] - beta1[len(beta1) - 1, :] # 计算残差
    62. res = np.array(res)
    63. press_i.append(np.sum(res ** 2)) # 残差平方和,并存储
    64. press.append(np.sum(press_i)) # 预测误差平方和
    65. # -----交叉有效性检验,判断主成分是否满足条件
    66. Q_h2[0] = 1
    67. if i > 0:
    68. Q_h2[i] = 1 - press[i] / ss[i - 1]
    69. if Q_h2[i] < 0.0975:
    70. print('提出的成分个数 r=', i + 1)
    71. break
    72. # -----根据主成分t计算回归方程的系数
    73. beta_Y_t = np.linalg.pinv(t) @ f0 # 求Y*关于t的回归系数
    74. beta_Y_X = w_star @ beta_Y_t # 求Y*关于X*的回归系数
    75. mu_x = mu[0:n] # 提取自变量的均值
    76. mu_y = mu[n:n + m] # 提取因变量的均值
    77. sig_x = sig[0:n] # 提取自变量的标准差
    78. sig_y = sig[n:n + m] # 提取因变量的标准差
    79. ch0 = mu_y - mu_x[:, np.newaxis].T / sig_x[:, np.newaxis].T @ beta_Y_X * sig_y[:, np.newaxis].T # 算原始数据回归方程的常数项
    80. beta_target = np.empty((n, 0)) # 回归方程的系数矩阵初始化
    81. for i in range(m):
    82. a = beta_Y_X[:, i][:, np.newaxis] / sig_x[:, np.newaxis] * sig_y[i] # 计算原始数据回归方程的系数
    83. beta_target = np.append(beta_target, a, axis=1)
    84. target = np.concatenate([ch0, beta_target], axis=0) # 回归方程的系数,每一列是一个方程,每一列的第一个数是常数项
    85. print(target)

    程序结果为:

    提出的成分个数 r= 2
    [[ 4.70197312e+01  6.12567103e+02  1.83984900e+02]
     [-1.66507920e-02 -3.50879763e-01 -1.25347701e-01]
     [-8.23702347e-01 -1.02476743e+01 -2.49692623e+00]
     [-9.69133092e-02 -7.41217596e-01 -5.18112643e-02]]

     3、时间序列分析

    ①自适应滤波法

    1. %% 自适应滤波法
    2. clc, clear
    3. yt = 0.1:0.1:1;
    4. m = length(yt);
    5. k = 0.9; % 学习常数
    6. N = 2; % 权数个数
    7. Terr = 10000;
    8. w = ones(1,N) / N;
    9. while abs(Terr)>0.00001
    10. Terr = [];
    11. for j = N+1:m-1
    12. yhat(j) = w*yt(j-1:-1:j-N)';
    13. err = yt(j) - yhat(j);
    14. Terr = [Terr, abs(Terr)];
    15. w = w + 2*k*err*yt(j-1:-1:j-N);
    16. end
    17. Terr = max(Terr);
    18. end

    ②指数平滑法

    二次指数平滑法

    1. %% 二次指数平滑法
    2. clc, clear
    3. yt = [676 825 774 716 940 1159 1384 1524 1668 1688 1958 2031 2234 2566 2820 3006 3093 3277 3514 3770 4107]';
    4. n = length(yt);
    5. alpha = 0.3;
    6. st1(1) = yt(1); st2(1) = yt(1); % 初始值选实际初始值
    7. for i = 2:n
    8. st1(i) = alpha*yt(i) + (1-alpha)*st1(i-1);
    9. st2(i) = alpha*st1(i) + (1-alpha)*st2(i-1);
    10. end
    11. at = 2*st1 - st2;
    12. bt = alpha/(1-alpha)*(st1-st2);
    13. yhat = at + bt;
    14. y_predict = alpha*yt(n) + (1-alpha)*yhat(n);
    15. x = 1965:1985;
    16. plot(x, yt, '*', x, yhat);

    ③移动平均法

    1、加权平均移动法

    1. %% 加权移动平均法
    2. clc ,clear
    3. y = [6.35 6.20 6.22 6.66 7.15 7.89 8.72 8.94 9.28 9.80];
    4. w = [1/6; 2/6; 3/6]; % 各加权项的权重
    5. m = length(y);
    6. n = 3; % 移动平均的项数
    7. for i = 1:m-n+1
    8. yhat(i) = y(i:i+n-1)*w;
    9. end
    10. err = abs(y(n+1:m)-yhat(1:end-1))./y(n+1:m); % 计算相对误差(最开始的三项数据是没有相对误差的)
    11. T_err = 1 - sum(yhat(1:end-1))/sum(y(n+1:m)); % 计算总的相对误差
    12. y_predict = yhat(end) / (1-T_err); % 用总的相对误差修正预测值

    2、趋势移动平均法

    1. %% 趋势移动平均法
    2. % 在一次移动平均的基础上再做一次移动平均
    3. clc, clear
    4. y = [676 825 774 716 940 1159 1384 1524 1668 1688 1958 2031 2234 2566 2820 3006 3093 3277 3514 3770 4107];
    5. % 第一次移动平均
    6. m1 = length(y);
    7. n = 6; % 移动平均的项数
    8. for i = 1:m1-n+1
    9. yhat1(i) = sum(y(i:i+n-1))/n;
    10. end
    11. % 第二次移动平均
    12. m2 = length(yhat1);
    13. for i = 1:m2-n+1
    14. yhat2(i) = sum(yhat1(i:i+n-1))/n;
    15. end
    16. plot(1:m1, y, '*');
    17. % 计算平滑系数
    18. at = 2*yhat1(end) - yhat2(end);
    19. bt = 2*(yhat1(end)-yhat2(end)) / (n-1);
    20. % 预测后两个数值
    21. y_predict1 = at + bt;
    22. y_predict2 = at + 2*bt;

    3、简单移动平均法

    1. %% 简单移动平均法
    2. clc, clear
    3. y = [533.8 574.6 606.9 649.8 705.1 772.0 816.4 892.7 963.9 1015.1 1102.7];
    4. m = length(y);
    5. n = [4,5]; % n位移动平均的项数,选取45分别做一次运算
    6. for i = 1:length(n)
    7. for j = 1:m-n(i)+1
    8. yhat{i}(j) = sum(y(j:j+n(i)-1))/n(i);
    9. end
    10. y_predict(i) = yhat{i}(end); % 预测值
    11. s(i) = sqrt(mean( (y(n(i)+1:m) - yhat{i}(1:end-1)).^2 )); % 预测的标准误差
    12. end

    4、微分方程

    参考之前的文章:《MATLAB--微分方程》

    http://t.csdn.cn/ub34W

    二、长期预测

    1、神经网络预测

    参考之前的文章:

    http://t.csdn.cn/zw7ou


     

    以及参考深度学习专栏

    搭建神经网络进行一些模型的预测

    2、logistic模型

    ①模型介绍

    • 优点:计算代价不高,易于理解和实现。
    • 缺点:容易欠拟合,分类精度可能不高。
    • 使用数据类型:数值型和标称型数据。

    ②模型分析及代码

    荷兰生物学家 Verhaust 于 1838 年提出了以昆虫数量为基础的 Logistic 人口增长模型。这个模型假设增长率 r 是人口的函数,它随着 x 的增加而减少。由 r(x)的表达式可知,当 x=𝑥𝑚时, r=0。其中, 𝑥𝑚表示自然资源条件能容纳的最大人口数。 Logistic 模型综合考虑了环境等因素对人口增长产生的影响,因此也是一种被广泛应用的效果较好的模型。\n\n设时刻 t 时人口为 x(t),环境允许的最大人口数量为𝑥𝑚,人口净增长率随人口数量的增加而线性减少,即:

     体现了环境和资源对人口增长的阻滞作用,显然 x 越大,该因子越大,说明随着人口数量 x 的增加,阻滞作用越来越大。\n\n由此建立 Logistic 人口模型的微分方程:

    根据附件给出的 1997-2016 年每年的人口数量,运用 MATLAB 进行编程处理,可以比较 1997-2016 年实际人口与拟合线的拟合程度,由图 1.1 可见, Logistic 模型对过去 20 年人口数量拟合效果较好。通过拟合,得到了待求参数 r=0.0394, 𝑥𝑚=1.5501× 105。因为以 1997 年为起算点,所以𝑥0=1.23626× 105。 则 Logistic 人口模型函数方程为:

     

     

     由 MATLAB 拟合出来的 Logistic 函数进行 2018-2068 年的人口预测, 人口数缓慢上升,经过 50 年人口增长了 13120 万人, Logistic 模型也求得在资源和环境的限制下, 无论过了多少年的人口不会超过 155010万人。在考虑资源和环境的阻滞作用随着人口增多也会越来越大的情况下,上述结果与实际是相符合的。

    1. clc,clear,
    2. t=1997:1:2016;
    3. x=[123626 124761 125786 126743 127627 128453 129227 129988 130756 131448
    4. 132129 132802 133450 134091 134735 135404 136072 136782 137462 138271];
    5. x1=[123626 124761 125786 126743 127627 128453 129227 129988 130756 131448
    6. 132129 132802 133450 134091 134735 135404 136072 136782 137462];
    7. x2=[124761 125786 126743 127627 128453 129227 129988 130756 131448 132129
    8. 132802 133450 134091 134735 135404 136072 136782 137462 138271];
    9. dx=(x2-x1)./x2;a=polyfit(x2,dx,1);
    10. r=a(2),xm=-r/a(1)
    11. x0=123626;
    12. f=inline('xm./(1+(xm/x0-1)*exp(-r*(t-1997)))','t','xm','r','x0');
    13. scatter(t,x,'k');
    14. hold on;
    15. plot(t,f(t,xm,r,x0),'r-');
    16. xlabel('年份');
    17. ylabel('人口数(万) ');
    18. title('1997-2016 年实际人口与拟合值比较');
    19. grid on;
    20. x2017=f(2017,xm,r,x0),x2018=f(2018,xm,r,x0),x2023=f(2023,xm,r,x0)
    21. x2028=f(2028,xm,r,x0),x2033=f(2033,xm,r,x0),x2038=f(2038,xm,r,x0)
    22. x2043=f(2043,xm,r,x0),x2048=f(2048,xm,r,x0),x2053=f(2053,xm,r,x0)
    23. x2058=f(2058,xm,r,x0),x2063=f(2063,xm,r,x0),x2068=f(2068,xm,r,x0)
    24. z=2018:5:2068;y=[x2018,x2023,x2028,x2033,x2038,x2043,x2048,x2053,x2058,x206
    25. 3,x2068];
    26. hold on;
    27. plot(z,y,'--r');xlabel('年份');ylabel('人口数(万) ');title('曲线拟合以及人口
    28. 预测曲线')
    29. grid on;
    1. clc,clear
    2. num=xlsread('C:\Users\\Desktop\Book1.xlsx','sheet1','B2:B21');
    3. x1=num;
    4. for i=1:18
    5. m1(i+2)=(x1(i)+x1(i+1)+x1(i+2))/3;
    6. end
    7. for i=3:18
    8. m2(i+2)=(m1(i)+m1(i+1)+m1(i+2))/3;
    9. end
    10. m1=m1';
    11. m2=m2';
    12. a=2*m1-m2;
    13. b=m1-m2;
    14. for i=1:52
    15. y(:,i)=a(20,1)+b(20,1)*i;
    16. end
    17. y
    18. figure;
    19. x=1:1:52;
    20. plot(x,y,'--r');
    21. xlabel('年份');
    22. ylabel('预测值');
    23. legend('预测值');
    24. title('2018-2068 年人口预测')

    文章应用案例来自:http://t.csdn.cn/9cmBB

  • 相关阅读:
    2019 国际大学生程序设计竞赛(ICPC)亚洲区域赛(银川) 7题
    <九> objectARX开发:读写Excel、json与txt格式文件
    计算有向图点的入度与出度
    bean的作用域、bean的生命周期、bean的后置处理器
    iOS程序内语言切换使用小结
    理想汽车 x JuiceFS:从 Hadoop 到云原生的演进与思考
    (附源码)ssm高校学生档案信息管理系统 毕业设计 010936
    卡尔曼滤波讲解
    总结sql用法及基础语法 第一章 三范式
    Spark底层原理:案例解析(第34天)
  • 原文地址:https://blog.csdn.net/m0_58585940/article/details/127720448