Matlab实现交通分布预测方法 —— 增长系数法 | 平均增长率法、底特律法、福莱特法

置顶 Douglassssssss

已于 2023-10-05 19:14:35 修改

阅读量3.8k

点赞数 18

分类专栏： Matlab 文章标签： matlab 算法交通物流笔记交通工程交通分布预测增长系数法

于 2023-05-02 11:28:39 首次发布

本文链接：https://blog.csdn.net/Douglassssssss/article/details/130459061

版权

Matlab 专栏收录该内容

2 篇文章 1 订阅

订阅专栏

前言

这是我的第一篇博客，其实我早就想写了，但每次都只是三分热情，一打开电脑就不知如何下手了，然后热情就没了。今天是五一假期第二天，还有两天假期，又不想看高数，又不想看专业课，但是又不想浪费时间，于是重拾了写博客的想法。

原来我以为可以用本地文件直接上传然后渲染的，现在看来还是得自己敲，不过还好CSDN是可以支持Md语法，也还算方便。

在用增长率法计算未来OD分布时，那几个方法用手算好繁琐，老师也讲了只是理解过程，最好还是得编程实现，我当然不放过这个机会。今天就来说一说如何用Matlab实现增长系数法预测未来OD交通量。

一、基本概念

交通分布预测：英文是Trip Distribution，是交通规划四阶段预测模型的第二步，是把交通的发生与吸引量预测获得的各小区的出行量转换成小区之间的空间OD量，即OD矩阵。

其实我现在看还是有点懵，不过看一下书上实际算例就能大致明白。

约束条件：所有小区的发生量之和等于总发生量，所有小区的吸引量之和等于总吸引量，总发生量等于总吸引量等于总交通生成量。

其实就是满足OD矩阵每行每列相加数字要对的上。

增长系数法：假定将来OD交通量的分布形式和现有的OD表的分布形式相同，在此假定的基础上预测对象区域目标年的OD交通量，常用的方法有常增长系数法、平均增长系数法、底特律（Detroit）法、福莱特（Fratar）法、佛尼斯法等。
发生、吸引增长率：未来发生（吸引）量和现状发生（吸引）量的比值，未来的是定值，现状的会随每一次迭代发生变化。
收敛判别：每一次迭代都会产生新的发生和吸引增长率，如果这两个增长率都在1附近浮动不大于3%（其他数也行，看精度要求），那么就可以停止迭代了。

二、实现过程

我们先来看看流程图。
算法流程图
图中m表示的是收敛需要的迭代次数，增长系数可以包含各种方法，因为不管是常增长还是福莱特法，都是计算一个系数乘上原来的OD矩阵，这也是为什么适合进行编程实现的原因。

我这里先附上代码，已经封装成了一个Matlab的function文件。

function [ res,o,d,m,T,Fo,Fd] = func_growthFactorMethod( X,U,V )
%FUNC_GROWTHFACTORMETHOD 增长系数法进行交通分布预测
%  res 为最终预测的od矩阵
%  o,d分别为未来各小区的总吸引量、发生量
%  m为迭代次数，T为未来出行生成总量
%  Fo,Fd分别为最终发生、吸引增长系数值
%  X 为现状OD矩阵 
%  U为未来各小区发生量，V为未来各小区吸引量，均为行向量。
n=size(X,1);     % 小区个数
O = sum(X,2)';    % 现状小区发生量
D = sum(X,1);    % 现状小区吸引量
% 计算发生增长系数和吸引增长系数
F_o = U ./ O;   F_d = V ./ D;
alpha = 0.03;    % 误差常数
flag = 0;        % 记录是否收敛
m = 0;     % 记录迭代次数
% 进行近似
% 选择近似方法
Choice = input('请选择方法（1-平均增长系数，2-底特律法，3-福莱特法）：  ');
if Choice == 1
   while true     % 开始循环 
       m = m+1;
       % 计算系数
       f = zeros(n);
       for i = 1:n
          for j = 1:n
             f(i,j)=0.5*(F_o(i)+F_d(j)); 
          end
       end
       % 计算新od表
       X = f .* X;
       % 更新现状O和D
       O = sum(X,2)';    
       D = sum(X,1);
       % 计算新的发生和吸引系数
       F_o = U ./ O;   F_d = V ./ D;
       % 判断收敛
       if abs([F_o,F_d]-1) < alpha
           flag = 1;
           res = X;
           o = sum(res,2)';
           d = sum(res,1);
           T = sum(o);
           Fo = F_o;
           Fd = F_d;
           break
       end
   end
end
if Choice == 2
   while true     % 开始循环 
       m = m+1;
       % 计算系数
       G = sum(O)/sum(U);  % 生成量增长率的倒数
       f = zeros(n);
       for i = 1:n
          for j = 1:n
             f(i,j)=G*F_o(i)*F_d(j); 
          end
       end
       % 计算新od表
       X = f .* X;
       % 更新现状O和D
       O = sum(X,2)';    
       D = sum(X,1);
       % 计算新的发生和吸引系数
       F_o = U ./ O;   F_d = V ./ D;
       G = sum(O)/sum(U);
       % 判断收敛
       if abs([F_o,F_d]-1) < alpha
           flag = 1;
           res = X;
           o = sum(res,2)';
           d = sum(res,1);
           T = sum(o);
           Fo = F_o;
           Fd = F_d;
           break
       end
   end
end
if Choice == 3
   while true     % 开始循环 
       m = m+1;
       % 计算系数
       f = zeros(n);
       for i = 1:n
          for j = 1:n
             Li = O(j)/sum((X(:,j)' .* F_d));
             Lj = D(i)/sum((X(i,:) .* F_o));
             f(i,j)=0.5*(Li+Lj)*F_o(i)*F_d(j); 
          end
       end
       % 计算新od表
       X = f .* X;
       % 更新现状O和D
       O = sum(X,2)';    
       D = sum(X,1);
       % 计算新的发生和吸引系数
       F_o = U ./ O;   F_d = V ./ D;
       % 判断收敛
       if abs([F_o,F_d]-1) < alpha
           flag = 1;
           res = X;
           o = sum(res,2)';
           d = sum(res,1);
           T = sum(o);
           Fo = F_o;
           Fd = F_d;
           break
       end
   end
end
end

我设置了一个用户交互式的输入，可以选择使用的方法（运行时需要在命令行输入对应的数字），佛尼斯法老师说不太用得到，我就没写。福莱特法比较难一些，需要注意的是位置系数应该是放在循环里的。

注意输入的U和V要写成行向量的形式（用逗号或空格隔开），如果写成了用分号隔开的列向量的形式也，输入的参数则应为其转置，即U’，V也是同理。下面给一个测试例子，也是书上的例子，方便大家复现。

原数据是这样的：

X = [17 7 4;
     7 38 6;
     4 5 17;];   % 现状od
U=[38.6 91.9 36]; % 未来各小区发生量
V=[39.3 90.3 36.9]; % 未来各小区吸引量
[ res,o,d,m,T] = func_growthFactorMethod(X, U, V);