2021华为杯E题思路+demo代码

2021 年中国研究生数学建模竞赛 E 题参考思路 交流群：912166339，非伸手党群 信号干扰下的超宽带（UWB）精确定位问题 一、背景 UWB（Ultra-Wideband）技术也被称之为“超宽带”，又称之为脉冲无线电技术。这是一 种无需任何载波，通过发送纳秒级脉冲而完成数据传输的短距离范围内无线通信技术，并且信 号传输过程中的功耗仅仅有几十µW。UWB 因其独有的特点，使其在军事、物联网等各个领域 都有着广阔的应用。其中，基于 UWB 的定位技术具备实时的室内外精确跟踪能力，定位精度 高，可达到厘米级甚至毫米级定位。UWB 在室内精确的定位将会对卫星导航起到一个极好的 补充作用，可在军事及民用领域有广泛应用，比如：电力、医疗、化工行业、隧道施工、危险 区域管控等。UWB 更多应用场景请参见[4—6]。 UWB 的定位技术有多种方法，本文仅考虑基于飞行时间（Time of Flight, TOF）的测距原 理，它是 UWB 定位法中最常见的定位方法之一。TOF 测距技术属于双向测距技术，其通过计 算信号在两个模块的飞行时间，再乘以光速求出两个模块之间的距离，这个距离肯定有不同程 度的误差，但其精度已经比较高。 在室内定位的应用中，UWB技术可以实现厘米级的定位精度（一般指2维平面定位）,并 具有良好的抗多径干扰和衰弱的性能以及具有较强的穿透能力。但由于室内环境复杂多变 UWB 通信信号极易受到遮挡，虽然UWB技术具有穿透能力，但仍然会产生误差，在较强干 扰时，数据会发生异常波动（通常是时间延时），基本无法完成室内定位，甚至会造成严重事 故。因此，信号干扰下的超宽带（UWB）精确定位问题成为亟待解决的问题。 二、问题描述 为解决信号干扰下的超宽带（UWB）精确定位问题，我们通过实际场景实测，采集到一 定数量的数据，即利用 UWB 的定位技术（TOF），采集到锚点（ anchor）与靶点（Tag）之间 的距离，希望通过数学建模（或算法）方法 ，无论信号是否干扰，都可以给出目标物（靶点） 的精确定位（ 3 维坐标）。 三、实验场景和数据采集 如图所示，在 5000mm5000mm3000mm 的测试环境中，分别在 4 个角落 A0，A1，A2， A3 放置 UWB 锚点（ anchor），锚点向所有方向发送信号。Tag 是 UWB 标签（靶点），即需 要定位的目标（只在测试环境范围内）。Tag 接收到 4 个 UWB 锚点（ anchor）的信号（无论 信号是否干扰，Tag 一般都可以接收到信号），利用 TOF 技术，分别解算出对应的 4 个距离数 据。 实验在实验场景 1 中采集了 Tag 在 324 个不同位置，在信号无干扰和信号干扰下的 UWB 数据，即每个位置各测试（采集）2 次，一次信号无干扰，另一次信号有干扰（锚点与靶点间 有遮挡），注意：每次采集数据时，由于 Tag 在同一位置会停留一会儿时间，而锚点与 Tag 之 间每 0.2—0.3 秒之间就会发送、接收信号一次，所以在同一位置点，UWB 会采集到多组数据 （多组数据都代表同一位置的信息），组数的多少视 Tag 在同一位置的时间而定，停留的时间 越长，组数就越多。数据见文件夹“附件 1：UWB 数据集”。 图 1 实测环境示意图 实验场景 1： 靶点（Tag）范围：5000mm5000mm3000mm 锚点（ anchor）位置（单位：mm）： A0（ 0，0，1300）、 A1（ 5000，0，1700）、 A2（ 0，5000，1700）、A3（ 5000，5000，1300） 四、数据文件说明 （ 1）UWB 数据集 “附件 1：UWB 数据集”有 2 个文件夹和 1 个文件，1 个文件（Tag 坐标信息.txt）存放 324 个不同位置的编号及 3 维坐标信息，2 个文件夹中 1 个存放信号无干扰下（正常）采集的 数据（各文件名为 x.正常.txt，x 表示对应的位置编号），另 1 个存放信号有干扰下（异常）采 集的数据（各文件名为 x.异常.txt，x 表示对应的位置编号）。 （ 2）数据文件 Tag 在每个位置都采集了 2 个数据文件（1 个正常，另 1 个异常），共有 648 个数据文件， 无论正常、异常数据，数据格式都一样，每个数据文件开头第 1 行为采集开始行，无实际意义， 接下来，每 4 行为一组，表示 UWB 采集的一组完整数据（一组数据表示一个样品），如： T:144235622:RR:0:0:950:950:118:1910 T:144235622:RR:0:1:2630:2630:118:1910 T:144235622:RR:0:2:5120:5120:118:1910 T:144235622:RR:0:3:5770:5770:118:1910 这 4 行数据的含义分别是: Tag 标识：时间戳：Range Report 的缩写：Tag ID：锚点 ID：该锚点的测距值(mm)：测 距值的校验值：数据序列号：数据编号（每个数据之间用“：”分隔）。实际上就是提供了 4 个锚点到靶点（Tag）的距离，即 A0 到靶点距离为：950mm A1 到靶点距离为：2630mm A2 到靶点距离为：5120mm A3 到靶点距离为：5770mm 每个数据文件都有多组数据，表示在同一位置连续时间内 UWB 自动采集到的多组数据。 五、完成任务研赛更多资料请加群：678535491，或关注公众号：老哥带你学术数模 试根据上述数据，完成如下任务： 任务 1：数据预处理（清洗） 无论是信号无干扰下采集数据，或信号干扰下采集数据，Tag 在同一坐标点上都采集多组 数据（见附件 1 中 648 个数据文件），请用某种方法把每个数据文件相应数值抓取出来，并转 换成二维表（矩阵）形式（ txt、Excel 或其他数据格式），每一行代表一组数据（即一个样品）， 然后对这些数据文件进行预处理（清洗），删除掉一些“无用”（异常、缺失、相同或相似） 的数据（样品）。经处理后，“正常数据”所有数据文件和“异常数据”所有数据文件最后各 保留多少组（多少个样品）数据，并重点列出以下 4 个数据文件，经处理后保留的数据（矩阵 形式）； “正常数据”文件夹中： 24.正常.txt、 109.正常.txt “异常数据”文件夹中： 1.异常.txt、 100.异常.txt 对于任务一而言，主要是进行数据的预处理工作；首先需要用 matlab 编写一套读取文件的代 码，参考如下： List =dir(‘C:\Users\DELL\Desktop\data*.txt’); %设置路径 for i=1:1:length(List) %挨个读入 txt 文件数据 A{i}=importdata([‘C:\Users\DELL\Desktop\data’,List(i).name]); End 然后对数据做相关的处理工作即可；可选择不同的方法将无用数据提取出来，其中异常数据可 以用简单统计分析、3σ原则或者箱线图来计算；对于缺失数据和相同及相似数据直接用程序 判断即可，找到上述异常数据后直接剔出。 3σ原则是建立在正态分布的等精度重复测量基础上而造成奇异数据的干扰或噪声难以满足正 态分布。研赛更多资料请加群：678535491，或关注公众号：老哥带你学术数模 正态分布：又叫做高斯分布。特征为中间高两边低左右对称。正态分布特性： 集中性：曲线的最高峰位于正中央，并且位置为均数所在的位置。 对称性：以均数所在的位置为中心呈左右对称，并且曲线两段无限趋近于横轴。 均匀变动性：正态分布曲线以均数所在的位置为中心均匀向左右两侧下降。 正态分布函数公式如下： 如果一组测量数据中的某个测量值的残余误差的绝对值 ν i>3σ，则该测量值为坏值，应剔除。 通常把等于 ±3σ的误差作为极限误差,对于正态分布的随机误差，落在 ±3σ以外的概率只 有 0.27%，它在有限次测量中发生的可能性很小，故存在 3σ准则。3σ准则是最常用也是最 简单的粗大误差判别准则，它一般应用于测量次数充分多( n ≥30)或当 n>10 做粗略判别时的 情况。 σ代表标准差,μ代表均值 样本数据服从正态分布的情况下： 数值分布在(μ-σ,μ+σ)中的概率为 0.6826 数值分布在(μ-2σ,μ+2σ)中的概率为 0.9544 数值分布在(μ-3σ,μ+3σ)中的概率为 0.9974 可以认为，Y 的取值几乎全部集中在（μ-3σ,μ+3σ)区间内，超出这个范围的可能性仅占不 到 0.3%。 任务 2: 定位模型 利用任务 1 处理后的数据，分别对“正常数据”和“异常数据”，设计合适的数学模型（或 算法），估计（或预测）出 Tag 的精确位置，并说明你所建立的定位模型（或算法）的有效性； 同时请利用你的定位模型（或算法）分别对附件 2 中提供的前 5 组（信号无干扰）数据和后 5 组（信号有干扰）数据进行精确定位（ 3 维坐标）； 注意：（1）定位模型必须体现实验场景信息； （ 2）请同时给出定位模型的 3 维（ x，y，z）精度、2 维（ x，y）精度以及 1 维的各自精度。 根据常识可知，若已知四个点的坐标，并且知道这四个点到另一个点的距离，此时能在空间上 确定该点的坐标；但题目中给了四个已知点，但存在大量的时间戳，因此存在多余解的问题， 对于多余解的问题则可以采用最小二乘的思想进行计算即可得到最终的点位坐标；而对于异常 数据，由于锚点和靶点间存在遮挡，因此此时计算的位置数据将存在偏差；此时可以先按照前 文正常数据计算的思路来时，将会发现计算得到的四组数据存在明显的偏差；因此可以首先将 可能认为是异常的点剔除（如统计所有样本数据将方差较大的值删除），然后利用其它正常点 来计算坐标即可。 另外在测量领域精度一定要用 3σ理念来表示，比较专业。 任务 3：不同场景应用 我们的训练数据仅采集于同一实验场景（实验场景 1），但定位模型应该能够在不同实际 场景上使用，我们希望你所建立的定位模型能够应用于不同场景。附件 3 中 10 组数据采集于 下面实验场景 2（前 5 组数据信号无干扰，后 5 组数据信号有干扰），请分别用上述建立的定 位模型，对这 10 组数据进行精确定位（ 3 维坐标）； 实验场景 2： 靶点（Tag）范围：5000mm3000mm3000mm 锚点（ anchor）位置（单位：mm）： A0（ 0，0，1200）、 A1（ 5000，0，1600）、 A2（ 0，3000，1600），A3（ 5000，3000，1200） 直接用任务 2 构建的模型进行计算即可 任务 4: 分类模型 上述定位模型是在已知信号有、无干扰的条件下建立的，但 UWB 在采集数据时并不知道 信号有无干扰，所以判断信号有无干扰是 UWB 精确定位问题的重点和难点。利用任务 1 处理 后的数据，建立数学模型(或算法)，以便区分哪些数据是在信号无干扰下采集的数据，哪些数 据是在信号干扰下采集的数据？并说明你所建立的分类模型(或算法)的有效性；同时请用你所 建立的分类模型（或算法）判断附件 4 中提供的 10 组数据（这 10 组数据同样采集于实验场景 1）是来自信号无干扰或信号干扰下采集的？ 这里是一道典型的分类问题，由于实现并不知道那些是有信号干扰的，那些是没有信号干扰的， 因此可以采用非监督聚类的思想进行分类工作，非监督聚类模型的主要有 K-means 聚类、高斯 混合聚类和层次聚类等，可以分别试试，类别写两类，这样最终就能将哪些数据是异常数据而 那些是正常数据进行划分；对于分类算法有效性判断可以利用附件 1：UWB 数据集进行验证， 任何分类算法都是又误差的，这时候需要把误差进行分析描述。如果误差较小则将其应用于附 件 4 中提供的 10 组数据即可。 K-means 算法代码 clear all; close all; clc; % 第一组数据 mu1=[0 0 ]; %均值(是需要生成的数据的均值) S1=[.1 0 ;0 .1]; %协方差(需要生成的数据的自相关矩阵（相关系数矩阵）) data1=mvnrnd(mu1,S1,100); %产生高斯分布数据 %第二组数据 mu2=[1.25 1.25 ]; S2=[.1 0 ;0 .1]; data2=mvnrnd(mu2,S2,100); % 第三组数据 mu3=[-1.25 1.25 ]; S3=[.1 0 ;0 .1]; data3=mvnrnd(mu3,S3,100); % 显示数据 plot(data1(:,1),data1(:,2),‘b+’); hold on;%不覆盖原图，要关闭则使用 hold off； plot(data2(:,1),data2(:,2),‘r+’); plot(data3(:,1),data3(:,2),‘g+’); grid on;%显示表格 % 三类数据合成一个不带标号的数据类 data=[data1;data2;data3]; N=3;%设置聚类数目 [m,n]=size(data);%表示矩阵 data 大小，m 行 n 列 pattern=zeros(m,n+1);%生成 0 矩阵 center=zeros(N,n);%初始化聚类中心 pattern(:,1:n)=data(:????; for x=1:N center(x,:)=data( randi(300,1)????;%第一次随机产生聚类中心 end while 1 %循环迭代每次的聚类簇； distence=zeros(1,N);%最小距离矩阵 num=zeros(1,N);%聚类簇数矩阵 new_center=zeros(N,n);%聚类中心矩阵 for x=1:m for y=1:N distence(y)=norm(data(x,:)-center(y,:));%计算到每个类的距离 end [~, temp]=min(distence);%求最小的距离 pattern(x,n+1)=temp;%划分所有对象点到最近的聚类中心；标记为 1,2,3； end k=0; for y=1:N for x=1:m if pattern(x,n+1)y new_center(y,:)=new_center(y,:)+pattern(x,1:n); num(y)=num(y)+1; end end new_center(y,:)=new_center(y,:)/num(y);%求均值，即新的聚类中心； if norm(new_center(y,:)-center(y,:))<0.1%检查集群中心是否已收敛。如果是则终止。 k=k+1; End end if kN break; else center=new_center; end end [m, n]=size(pattern); %最后显示聚类后的数据 figure; hold on; for i=1:m if pattern(i,n)1 plot(pattern(i,1),pattern(i,2),‘r*’); plot(center(1,1),center(1,2),‘ko’);%用小圆圈标记中心点； elseif pattern(i,n)2 plot(pattern(i,1),pattern(i,2),‘g*’); plot(center(2,1),center(2,2),‘ko’); elseif pattern(i,n)3 plot(pattern(i,1),pattern(i,2),‘b*’); plot(center(3,1),center(3,2),‘ko’); elseif pattern(i,n)4 plot(pattern(i,1),pattern(i,2),‘y*’); plot(center(4,1),center(4,2),‘ko’); else plot(pattern(i,1),pattern(i,2),‘m*’); plot(center(4,1),center(4,2),‘ko’); end end grid on; 高斯混合聚类代码 %% 初始化工作空间 clc clear all close all %% 载入数据 load fisheriris %% 二维数据 % 叶子长度和叶子宽度散点图(真实标记) figure; speciesNum = grp2idx(species); gscatter(meas(:,3),meas(:,4),speciesNum,[‘r’,‘g’,‘b’]) xlabel(‘叶子长度’) ylabel(‘叶子宽度’) title(‘真实标记’) set(gca,‘FontSize’,12) set(gca,‘FontWeight’,‘bold’) % 叶子长度和叶子宽度散点图（无标记） figure; scatter(meas(:,3),meas(:,4),150,’.’) xlabel(‘叶子长度’) ylabel(‘叶子宽度’) title(‘无标记’) set(gca,‘FontSize’,12) set(gca,‘FontWeight’,‘bold’) %% 高斯混合模型聚类 data = [meas(:,3), meas(:,4)]; % 手动初始条件 Mu=[0.25 1.5; 4.0 1.25; 5.5 2.0 ];%平均值设定 Sigma(:,:,1) = [1 1;1 2];%方差设定 Sigma(:,:,2) = [1 1;1 2]; Sigma(:,:,3) = [1 1;1 2]; Pcom=[1/3 1/3 1/3];%混合比例 S = struct(‘mu’,Mu,‘Sigma’,Sigma,‘ComponentPropotion’,Pcom); GMModel = fitgmdist(data,3,‘Start’,S); % 分类 T1 = cluster(GMModel,data); % 标号调整 cen=[mean(data(T11,:));… mean(data(T12,:));… mean(data(T13,:))]; dist=sum(cen.^2,2); [dump,sortind]=sort(dist,‘ascend’); newT1=zeros(size(T1)); for i =1:3 newT1(T1i)=find(sortindi); end % 叶子长度和叶子宽度散点图(真实标记:实心圆+kmeans 分类:圈) figure; gscatter(meas(:,3),meas(:,4),speciesNum,[‘r’,‘g’,‘b’]) hold on gscatter(data(:,1),data(:,2),newT1,[‘r’,‘g’,‘b’],‘o’,10) scatter(cen(:,1),cen(:,2),300,‘m*’) hold off xlabel(‘叶子长度’) ylabel(‘叶子宽度’) title(‘真实标记:实心圆+kmeans 分类:圈’) set(gca,‘FontSize’,12) set(gca,‘FontWeight’,‘bold’) %% 混淆矩阵 ConfusionMatrix T2ConfMat=confusionmat(speciesNum,newT1) error23=(speciesNum2)&(newT13); errDat23=data(error23,:) error32=(speciesNum3)&(newT1==2); errDat32=data(error32,:) %% 高斯模型等高线图 % 散点图 figure; gscatter(meas(:,3),meas(:,4),speciesNum,[‘r’,‘g’,‘b’]) hold on scatter(cen(:,1),cen(:,2),300,‘m*’) hold off xlabel(‘叶子长度’) ylabel(‘叶子宽度’) title(‘高斯模型等高线图’) set(gca,‘FontSize’,12) set(gca,‘FontWeight’,‘bold’) % 叠加等高线 haxis=gca; xlim = get(haxis,‘XLim’); ylim = get(haxis,‘YLim’); dinter=(max([xlim, ylim]) - min([xlim, ylim]))/100; [Grid1, Grid2] = meshgrid(xlim(1):dinter:xlim(2), ylim(1):dinter:ylim(2)); hold on GMMpdf=reshape(pdf(GMModel, [Grid1(???? Grid2(????]), size(Grid1,1), size(Grid2,2)); contour(Grid1, Grid2, GMMpdf, 30); hold off % 混合高斯模型曲面图 figure; surf(GMMpdf) xlabel(‘叶子长度’) ylabel(‘叶子宽度’) title(‘高斯混合模型曲面图’) view(-3,65) 层次聚类相关代码 %% 读取数据 clc;clear all,close all pre_data=xlsread(‘StockFinance.xlsx’,‘Sheet1’,‘E2:N4735’); %% 数据归一化 [rn,cn]=size(pre_data); step1_data=zeros(rn,cn); for k=1:cn %基于均值方差的离群点数据归一化 xm=mean(pre_data(:,k)); xs=std(pre_data(:,k)); for j=1:rn if (pre_data(j,k))>xm+2xs step1_data(j,k)=1; elseif (pre_data(j,k))xs step1_data(j,k)=0; else step1_data(j,k)=(pre_data(j,k)-(xm-2xs))/(4xs); end end end xlswrite(‘norm_data.xlsx’,step1_data); %% 层次聚类 numClust=3; dist_h=‘spearman’; link=‘weighted’; hidx=clusterdata(step1_data,‘maxclust’,numClust,‘distance’,dist_h,‘linkage’,link); 任务 5: 运动轨迹定位 运动轨迹定位是 UWB 重要应用之一，利用静态点的定位模型，加上靶点自身运动规律， 希望给出动态靶点的运动轨迹。附件 5 是对动态靶点采集的数据（一段时间内连续采集的多组 数据），请注意，在采集这些数据时，会随机出现信号干扰，请对这个运动轨迹进行精确定位， 最终画出这条运动轨迹图（数据采集来自实验场景 1）。 这个就是结合前面各小问的解决方案进行的集成，首先题目明确告知了会随机出现干扰，，也 就意味着在整个时间线内即存在正常数据也存在干扰数据，因此因该先用附件 4 中的聚类模型 将干扰数据进行筛选，然后利用任务 2 构建的基于干扰数据的定位模型进行点位坐标的确定， 再结合正常数据的坐标即可得到最终的运动轨迹信息。 参考文献 [1]Kok M, Hol J D, Schon T B. Indoor Positioning Using Ultrawideband and Inertial [2]Measurements[J]. IEEE Transactions on Vehicular Technology, 2015, 64(4): 1293—1303. [2] Arias-De-Reyna E. A Cooperative Localization Algorithm for UWB Indoor Sensor Networks[J]. Wireless Personal Communications, 2013, 72(1):85—99. [3]缪希仁,范建威等，基站异常情况下基于改进极限学习机的超宽带室内定位方法，传感技术 学报，2020，33（ 10）：1—10. [4]http://m.elecfans.com/article/1643155.html [5]http://m.elecfans.com/article/1116045.html [6]http://www.qianjia.com/zhike/html/2020-02/11_19601.html