匈牙利算法与卡尔曼滤波
匈牙利算法与卡尔曼滤波
- 匈牙利算法
- 指派问题概述
- 匈牙利算法流程
- 卡尔曼滤波
- 一维卡尔曼滤波
- 卡尔曼滤波器方程
- 总结
- 多维卡尔曼滤波
- 扩展卡尔曼滤波器（EKF）
- 扩展卡尔曼滤波器【EKF】雅可比矩阵
- UFK 无迹卡尔曼滤波
匈牙利算法
指派问题概述
首先，对匈牙利算法解决的问题进行概述：实际中，会遇到这样的问题，有n项不同的任务，需要n个人分别完成其中的1项，每个人完成任务的时间不一样。于是就有一个问题，如何分配任务使得花费时间最少。
通俗来讲，就是n*n矩阵中，选取n个元素，每行每列各有一个元素，使得和最小。
可以抽象成一个矩阵，如果是求和最小问题，那么这个矩阵就叫做花费矩阵（Cost Matrix）；如果要求的问题是使之和最大化，那么这个矩阵就叫做利益矩阵（Profit Matrix）.
匈牙利算法流程
匈牙利算法首先进行列规约，即每行减去此行最小元素，每一列减去该列最小元素，规约后每行每列中必有0元素出现。接下来进行试指派，也就是划最少的线覆盖矩阵中全部的0元素，如果试指派的独立0元素数等于方阵维度则算法结束，如果不等于则需要对矩阵进行调整，重复式指派和调整步骤直到满足算法结束条件。
import numpy as np
import collections
import timeclass Hungarian():def __init__(self, input_matrix=None, is_profit_matrix=False):"""输入为一个二维嵌套列表is_profit_matrix = False代表输入是消费矩阵（需要使消费最小化），反之则为利益矩阵（需要使利益最大化）"""if input_matrix is not None:# 保存输入my_matrix = np.array(input_matrix)self._input_matrix = np.array(input_matrix)self._maxColumn = my_matrix.shape[1]self._maxRow = my_matrix.shape[0]# 本算法必须作用于方阵，如果不为方阵则填充为方阵matrix_size = max(self._maxColumn, self._maxRow)pad_columns = matrix_size - self._maxRowpad_rows = matrix_size - self._maxColumnmy_matrix = np.pad(my_matrix, ((0, pad_columns), (0, pad_rows)), 'constant', constant_values=(0))# 如果需要，则转化为消费矩阵if is_profit_matrix:my_matrix = self.make_cost_matrix(my_matrix)self._cost_matrix = my_matrixself._size = len(my_matrix)self._shape = my_matrix.shape# 存放算法结果self._results = []self._totalPotential = 0else:self._cost_matrix = Nonedef make_cost_matrix(self, profit_matrix):'''利益矩阵转化为消费矩阵，输出为numpy矩阵'''# 消费矩阵 = 利益矩阵最大值组成的矩阵 - 利益矩阵matrix_shape = profit_matrix.shapeoffset_matrix = np.ones(matrix_shape, dtype=int) * profit_matrix.max()cost_matrix = offset_matrix - profit_matrixreturn cost_matrixdef get_results(self):"""获取算法结果"""return self._resultsdef calculate(self):"""实施匈牙利算法的函数"""result_matrix = self._cost_matrix.copy()# 步骤1：矩阵每一行减去本行的最小值for index, row in enumerate(result_matrix):result_matrix[index] -= row.min()# 步骤2：矩阵每一列减去本行的最小值for index, column in enumerate(result_matrix.T):result_matrix[:, index] -= column.min()# 步骤3：使用最少数量的划线覆盖矩阵中所有的0元素# 如果划线总数不等于矩阵的维度，需要进行矩阵调整并重复循环此步骤total_covered = 0while total_covered < self._size:time.sleep(1)# 使用最少数量的划线覆盖矩阵中所有的0元素同时记录划线数量cover_zeros = CoverZeros(result_matrix)single_zero_pos_list = cover_zeros.calculate()covered_rows = cover_zeros.get_covered_rows()covered_columns = cover_zeros.get_covered_columns()total_covered = len(covered_rows) + len(covered_columns)# 如果划线总数不等于矩阵的维度需要进行矩阵调整（需要使用未覆盖处的最小元素）if total_covered < self._size:result_matrix = self._adjust_matrix_by_min_uncovered_num(result_matrix, covered_rows, covered_columns)# 元组形式结果对存放到列表self._results = single_zero_pos_list# 计算总期望结果value = 0for row, column in single_zero_pos_list:value += self._input_matrix[row, column]self._totalPotential = valuedef get_total_potential(self):return self._totalPotentialdef _adjust_matrix_by_min_uncovered_num(self, result_matrix, covered_rows, covered_columns):"""计算未被覆盖元素中的最小值（m）,未被覆盖元素减去最小值m,行列划线交叉处加上最小值m"""adjusted_matrix = result_matrix# 计算未被覆盖元素中的最小值（m）elements = []for row_index, row in enumerate(result_matrix):if row_index not in covered_rows:for index, element in enumerate(row):if index not in covered_columns:elements.append(element)min_uncovered_num = min(elements)#print('min_uncovered_num:',min_uncovered_num)#未被覆盖元素减去最小值mfor row_index, row in enumerate(result_matrix):if row_index not in covered_rows:for index, element in enumerate(row):if index not in covered_columns:adjusted_matrix[row_index,index] -= min_uncovered_num#print('未被覆盖元素减去最小值m',adjusted_matrix)#行列划线交叉处加上最小值mfor row_ in covered_rows:for col_ in covered_columns:#print((row_,col_))adjusted_matrix[row_,col_] += min_uncovered_num#print('行列划线交叉处加上最小值m',adjusted_matrix)return adjusted_matrixclass CoverZeros():"""使用最少数量的划线覆盖矩阵中的所有零输入为numpy方阵"""def __init__(self, matrix):# 找到矩阵中零的位置（输出为同维度二值矩阵，0位置为true，非0位置为false）self._zero_locations = (matrix == 0)self._zero_locations_copy = self._zero_locations.copy()self._shape = matrix.shape# 存储划线盖住的行和列self._covered_rows = []self._covered_columns = []def get_covered_rows(self):"""返回覆盖行索引列表"""return self._covered_rowsdef get_covered_columns(self):"""返回覆盖列索引列表"""return self._covered_columnsdef row_scan(self,marked_zeros):'''扫描矩阵每一行，找到含0元素最少的行，对任意0元素标记（独立零元素），划去标记0元素（独立零元素）所在行和列存在的0元素'''min_row_zero_nums = [9999999,-1]for index, row in enumerate(self._zero_locations_copy):#index为行号row_zero_nums = collections.Counter(row)[True]if row_zero_nums < min_row_zero_nums[0] and row_zero_nums!=0:#找最少0元素的行min_row_zero_nums = [row_zero_nums,index]#最少0元素的行row_min = self._zero_locations_copy[min_row_zero_nums[1],:]#找到此行中任意一个0元素的索引位置即可row_indices, = np.where(row_min)#标记该0元素#print('row_min',row_min)marked_zeros.append((min_row_zero_nums[1],row_indices[0]))#划去该0元素所在行和列存在的0元素#因为被覆盖，所以把二值矩阵_zero_locations中相应的行列全部置为falseself._zero_locations_copy[:,row_indices[0]] = np.array([False for _ in range(self._shape[0])])self._zero_locations_copy[min_row_zero_nums[1],:] = np.array([False for _ in range(self._shape[0])])def calculate(self):'''进行计算'''#储存勾选的行和列ticked_row = []ticked_col = []marked_zeros = []#1、试指派并标记独立零元素while True:#print('_zero_locations_copy',self._zero_locations_copy)#循环直到所有零元素被处理（_zero_locations中没有true）if True not in self._zero_locations_copy:breakself.row_scan(marked_zeros)#2、无被标记0（独立零元素）的行打勾independent_zero_row_list = [pos[0] for pos in marked_zeros]ticked_row = list(set(range(self._shape[0])) - set(independent_zero_row_list))#重复3,4直到不能再打勾TICK_FLAG = Truewhile TICK_FLAG:#print('ticked_row:',ticked_row,' ticked_col:',ticked_col)TICK_FLAG = False#3、对打勾的行中所含0元素的列打勾for row in ticked_row:#找到此行row_array = self._zero_locations[row,:]#找到此行中0元素的索引位置for i in range(len(row_array)):if row_array[i] == True and i not in ticked_col:ticked_col.append(i)TICK_FLAG = True#4、对打勾的列中所含独立0元素的行打勾for row,col in marked_zeros:if col in ticked_col and row not in ticked_row:ticked_row.append(row)FLAG = True#对打勾的列和没有打勾的行画画线self._covered_rows = list(set(range(self._shape[0])) - set(ticked_row))self._covered_columns = ticked_colreturn marked_zeros
卡尔曼滤波
现代系统大多数都配备了多个传感器，这些传感器根据一系列测量值提供隐藏变量的估计。例如，GPS接收机提供位置和速度估计，其中位置和速度是隐藏变量，卫星信号到达差值时间是测量值。
一维卡尔曼滤波
系统状态 当前状态是预测算法的输入，下一个状态（下一时间间隔的目标参数）是算法输出。
动态模型 描述输入和输出之间的关系
测量噪声 测量中包含的误差
过程噪声 动态模型误差
https://zhuanlan.zhihu.com/p/80255855
卡尔曼滤波器方程
- 状态更新方程 State Update Equation
这里的factor成为**卡尔曼增益 （Kalman Gain）*用KnK_nKn表示。下标n表示卡尔曼增益随每次迭代而变化。
也就是说，当前状态等于卡尔曼滤波在n时刻的预测值＋系数（n时刻测量值-n时刻的预测值） - 状态外推方程（过渡方程或预测方程）
这个方程组将当前状态外推到下一个状态，即判断状态是否更新，或怎么更新。 - 卡尔曼增益方程
当测量不确定度非常大且估计不确定度非常小时，卡尔曼收益趋近于0.因此，我们估计值一个较大的权重，给测量值一个较小的权重。
卡尔曼增益>>通过给定的测量值改变我的估计。
卡尔曼增益越低，越依赖估计值；相反，越依赖测量值。 - 协方差更新方程
由于（1-Kn）<=1，我们可以很清楚的从方程中看到，随着滤波器的每次迭代，估计不确定行会变得越来越小。当测量不确定性增大，卡尔曼增益会变小，因此，估计不确定度的收敛会变慢。 - 协方差外推方程
利用动态模型方程进行。
总结
- 初始化时，当估计误差很大时，估计初值并不重要，因为这样会计算出较大的卡尔曼增益导致预测值与测量值非常接近。
- 如果我们使用了一个恒定系统动力模型的一维卡尔曼滤波器来测量温度变化的液体温度。我们在卡尔曼滤波器估计中观察了滞后误差。这个滞后误差是由错误的系统动力模型与错误的过程模型导致的。（ps:滞后误差可以通过合理的定义系统动力模型与过程模型得到解决）
- 一个较大的过程误差，能够维持协方差不会一致缩小，使得卡尔曼增益维持一个较大的水平，从而更加贴近观测值。
多维卡尔曼滤波
卡尔曼滤波的五个方程分别为：状态方程、协方差方程，状态更新方程、协方差更新方程、卡尔曼增益方程。
B：表示控制矩阵
Q：过程噪声
F：转移矩阵
z: 观测值
R： 观测方差
HPH^T: 估计方差
扩展卡尔曼滤波器（EKF）
处理非线性问题
有两件事要考虑：
- 在不知道它的基本函数的情况下，如何从实际信号中计算一阶导数；
- 如何将我们的单值非线性状态/观测模型推广到我们一直在考虑的多值系统中；
如果一个信号（如传感器值z_k）是另一个信号的函数（如状态x_k），我们可以使用第一个信号的连续差除以第二个信号的连续差：
扩展卡尔曼滤波器【EKF】雅可比矩阵
雅可比矩阵
如何将我们的单值非线性状态/观测模型推广到多值系统，这将有助于我们回忆线性模型传感器组件的方程：
zk=Cxkz_k=Cx_kzk=Cxk
对于一个有两个状态值与三个传感器的系统，我们可以重写系统观测模型：
对于非线性模型，同样会有一个矩阵，其行数等于传感器数，列数等于状态数；但是，该矩阵将包含传感器值相对于该状态值的一阶导数的当前值。数学家将这种导数称为偏导数，并将这些导数的矩阵称为雅可比矩阵。
例如一项调查，其中一组人被要求在一个量表上对两种不同的产品进行评分。每个产品的总分将是该产品所有用户评分的平均值。为了了解一个人如何影响单个产品的整体评级，我们将查看此人对该产品的评级。每个人都像是一个偏导数，而这个人就像是雅可比矩阵。用传感器替换人，用状态替换问题，你就能理解扩展卡尔曼滤波器的传感器模型。
我们的线性模型：
xk=Axk−1+wkx_k=Ax_{k-1}+w_kxk=Axk−1+wk
变成：
xk=f(xk−1)+wkx_k=f(x_{k-1})+w_kxk=f(xk−1)+wk
在非线性的情况下，我们可以使用扩展卡尔曼滤波器（EKF），它把非线性函数在当前估计状态的平均值附近进行线性化。
在每个微小时刻执行线性化，然后将得到的雅可比矩阵用于预测和更新卡尔曼滤波器的算法状态。
当系统是非线性的，并且可以通过线性化很好的近似时，那么扩展卡尔曼滤波器（EKF）是状态估计的一个很好的选择。
但是，扩展卡尔曼滤波器有以下缺点：
- 由于复杂的导数，可能难以解析计算雅可比矩阵
- 以数值方式计算雅可比矩阵，可能需要很高的计算成本
- 扩展卡尔曼滤波器不适用于不连续的系统，因为系统不可微分时，雅可比矩阵不存在
- 对于高度非线性化的系统，效果并不好。
UFK 无迹卡尔曼滤波
逼近概率分布要比逼近任意的非线性函数要容易得多。
匈牙利算法与卡尔曼滤波相关推荐
- 目标跟踪--卡尔曼滤波 与 匈牙利算法
目前主流的目标跟踪算法都是基于Tracking-by-Detecton策略,即基于目标检测的结果来进行目标跟踪. 跟踪结果中,每个bbox左上角的数字是用来标识某个人的唯一ID号.那么问题就来了,视频 ...
- 车流量检测实现：多目标追踪、卡尔曼滤波器、匈牙利算法、SORT/DeepSORT、yoloV3、虚拟线圈法、交并比IOU计算
日萌社 人工智能AI:Keras PyTorch MXNet TensorFlow PaddlePaddle 深度学习实战(不定时更新) CNN:RCNN.SPPNet.Fast RCNN.Faste ...
- 匈牙利匹配、KM算法、卡尔曼滤波、SORT/deep SORT
匈牙利匹配.KM算法 带你入门多目标跟踪(三)匈牙利算法&KM算法 - 知乎 卡尔曼滤波: 图像处理之目标跟踪(一)之卡尔曼kalman滤波跟踪(主要为知识梳理)(转载)_coming_is_ ...
- 多目标跟踪卡尔曼滤波和匈牙利算法
多目标跟踪关联匹配算法(匈牙利算法和KM算法原理讲解和代码实现) 目录 多目标跟踪关联匹配算法 多目标跟踪关联匹配算法(匈牙利算法和KM算法原理讲解和代码实现) 0.多目标跟踪算法流程 1.卡尔曼滤波 ...
- 目标跟踪中的卡尔曼滤波和匈牙利算法解读。
先解读Sort算法:Simple online and realtime tracking 论文地址 https://arxiv.org/abs/1602.00763 代码地址 https://git ...
- 匈牙利算法和卡尔曼滤波器
1.匈牙利算法 **分配问题(Assignment Problem):**假设有N个人和N个任务,每个任务可以任意分配给不同的人,已知每个人完成每个任务要花费的代价不尽相同,那么如何分配可以使得总的代 ...
- 匈牙利算法 求二分图最大匹配
匈牙利算法 1. 二分图 二分图: 又称作二部图,是图论中一种特殊模型.设G=(V,E)是一个无向图,如果顶点V可分割为两个互不相交的子集(A,B),并且图中每条边所关联的两个顶点 i 和 j 分别属 ...
- 二分图匹配匈牙利算法DFS实现
1 /*==================================================*\ 2 | 二分图匹配(匈牙利算法DFS 实现) 3 | INIT: g[][]邻接矩阵; ...
- 解题报告：luogu P2423 [HEOI2012]朋友圈【最大团转最大点独立集（匈牙利算法+时间戳优化）】
图的最大团:"任意两点之间都有一条边相连"的子图被称为无向图的团,点数最多的团为图的最大团 朋友圈中任意两个点之间都有关系,既是图中的团. 答案就是图中的最大团. 我们如果把B国的 ...
最新文章
- 微信小程序——云服务环境的配置
- 判断101-200之间有多少个素数，并输出所有素数。
- 一个用户的上级部门的上级部门对用户也有修改权限，怎么判断？
- skills --札记
- 零基础学python用哪本书好-零基础学习python推荐几本书？
- uc浏览器邀请码_UC密保手机不能用？冬树教你如何一招申诉成功！
- 杭电2006~2009计算机学院笔试真题详解
- html禁止页面动画,Animate.css 超强CSS3动画库，三行代码搞定H5页面动画特效！
- java 虚拟机内存修改_Java虚拟机内存参数设置
- 步进驱动器简单接线说明书
- python实现r树存储地理位置_空间索引R树研究_回顾与展望_张明波
- 最后两星期，怎么过6级？（最快攻略）
- 利用INFOPATH2007VS2005开发MOSS工作流详解 --收藏
- three.js 渲染调优，如何提升3d场景更逼真的渲染效果
- 前端面试宝典React篇03 如何避免生命周期中的坑？
- Discussion 2
- 拼多多求变 200 天：撒钱百亿元，江湖人称拼爹爹？
- 【181018】基于MFC文档方式制作的飞碟射击游戏
- 15.1 DIB 文件格式
- 查询自己电脑的IP地址