【项目介绍】基于机器学习的低空小、微无人机识别技术

1.项目介绍

对于现代雷达探测系统而言，无人机和飞鸟同属于低空小、微特征的一类典型目标，而面对比较复杂的环境，如何有效区分两者类型并完成识别是当下急迫且重要的难题。常规方法是从目标的微动特征差异进行区分，但由于两者回波微弱，很难通过时频分析方法提取目标特征。针对无人机与鸟类轨迹的特性差异，我们设计了多维特征提取方法，包括轨迹角度变化、航向角振荡、速度分布等物理量，为分类模型提供了充分的信息支持。我们采用了多种主流机器学习算法的组合，通过Stacking集成学习方法，有效提升了模型的预测能力和泛化性能。本项目聚焦无人机与鸟类航迹识别的挑战，旨在通过深度分析两者在运动特性上的差异，开发一个基于特征提取和机器学习模型的低空小、微无人机识别技术。无人机和鸟类的飞行轨迹具有显著不同的动态特性，如何准确区分二者是智能监控、边境安防、航空安全、环境与生态监测等场景中的重要任务。为此，本项目采用了一套严谨的预处理与特征提取流程，并通过先进的集成学习方法提升模型的精度和泛化能力。首先对比两者在运动轨迹等信息上的差异性，进行特征分析，提出了时间相关的航向震荡频率与速度分布等特征量描述方法，利用雷达系统记录的航迹数据，提取两者的有效特征量；然后利用XGBoost算法对样本进行训练，并在获得最优模型参数后，通过测试样本进行测试，测试分类结果显示识别的准确率能够达93%，其结果表明了该方法的正确性也体现了该方法在识别中的稳定性和可靠性，具有较高的实际应用价值。

2.数据预处理

项目原始数据包含无人机和飞鸟的雷达数据信息，标签为1和0，但并不能确定1和0分别是哪一类。先将这两类的轨迹画出来，观察轨迹图，因为无人机和飞鸟的轨迹跨度较大，所以先将其统一轨迹范围后再画图，从视觉上先感受无人机合和飞鸟轨迹的区别。
由于初始数据含有较多重复数据和异常数据，所以我们需要先对初始数据进行数据清洗,数据一般在时间上都是间隔2s的，但是在一组无人机/飞鸟的数据内，可能会出现时间上的不连续，所以我们预处理代码先把这些时间间隔出现突变的数据，在标签列标记上3，在以后读取的时候，不能将含有3的这一组数据视为1组数据。下面的代码是去除数据中的重复文件

% 此文件检测数据文件中的重复值，并删除，保留最初一行，保持绝对行不变
clear;
clc;

% 读取Excel数据时忽略首行
data = readtable("data.xlsx", 'ReadVariableNames', false, 'Range', 'A:G');


% 提取第7列数据
column7_data = data{:, 7};

% 初始化一个空的单元数组，大小和 column7_data 相同


% 遍历每个 cell 并进行处理
% for i = 1:length(column7_data)
%     if isempty(column7_data{i})
%         % 如果 cell 为空，保留空值（可以选择 NaN 或其他占位符）
%         output_data{i} = NaN;  % 对于数值数据，用 NaN 作为占位符
%     elseif isnumeric(column7_data{i})
%         % 如果是数值，直接保留
%         output_data{i} = column7_data{i};
%     elseif ischar(column7_data{i})
%         % 如果是字符串，转换为数值占位符，或选择其他处理方法
%         output_data{i} = output_data{i};  % 或者可以选择不同的占位符，如 0 或空字符串
%     else
%         % 处理其他类型的数据
%         output_data{i} = NaN;  % 可根据需要替换
%     end
% end
% 
% 
% % 最后将 cell 转换为数值矩阵
% output_data = cell2mat(output_data);
% output_data = array2table(output_data);
% 
% 
% 
% data(:,7) = [];
% data(:,7) = output_data;

if iscell(column7_data)
    column7_data = str2double(column7_data);  % 转换为数值
end
column7_data(end+1,:) = 0;
output_data = cell(size(column7_data));
% 提取数据
azimuth = data{:, 1};    % 方位角
range = data{:, 2};      % 距离
relative_height = data{:, 3};  % 相对高度
labels = data{:, 7};     % 标签列
RCS = 1000 * data{:, 6}; % RCS
vel_jingxiang = data{:, 4}; % 径向速度
time = data{:, 5};

% 检查标签列是否为cell类型，如果是，则将其转换为合适的数值格式
if iscell(labels)
    labels = cellfun(@str2double, labels);  % 使用 cellfun 转换为数值
end

% 无人机记录起始点和终止点索引
start_indices = find(~isnan(column7_data));  % 找到标签为0的起始点

% 初始化保存所有data_seg的表格
all_data_segs = table();  % 保留为table格式

for i = 1:length(start_indices) - 1
    current_start = start_indices(i);
    current_end = start_indices(i + 1) - 1;

    % 提取该时间段内的数据
    data_seg = data(current_start:current_end, :);
    data_without_labels = data_seg(:, 1:6);

    % 删除重复行
    [~, unique_idx] = unique(data_without_labels, 'rows', 'stable');
    duplicate_rows = setdiff(1:size(data_without_labels, 1), unique_idx);
    num_duplicates = length(duplicate_rows);

    % 删除重复行
    data_seg(duplicate_rows, :) = [];

    % 生成具有相同变量名的空行（保持表格类型一致）
    empty_rows = table(zeros(length(duplicate_rows),1)*NaN, ...
        zeros(length(duplicate_rows),1)*NaN, ...
        zeros(length(duplicate_rows),1)*NaN,...
        zeros(length(duplicate_rows),1)*NaN, ...
        zeros(length(duplicate_rows),1)*NaN, ...
        zeros(length(duplicate_rows),1)*NaN, ...
        zeros(length(duplicate_rows),1)*NaN);
    

    if num_duplicates > 0
    % 只有当 num_duplicates > 0 时才执行这行代码
    data_seg = [data_seg; empty_rows];
    end

    
    % 将 data_seg 追加到 all_data_segs
    all_data_segs = [all_data_segs; data_seg];
    %disp("执行了%d\n",i)
end

% disp(all_data_segs)
% 保存所有 data_seg 到 Excel 文件
writetable(all_data_segs, 'data_test_chongfu.xlsx');

经过去重操作后，给时间超出阈值的数据标记为3

clear;

% 读取Excel数据
data = readtable('data_test_chongfu.xlsx');

% 提取数据列
azimuth = data{:, 1};    % 方位角
range = data{:, 2};      % 斜距
relative_height = data{:, 3};  % 相对高度
v_h = data{:, 4};        % 径向速率
time = data{:, 5};       % 时间列
RCS = data{:, 6};        % RCS
labels = data{:, 7};     % 标签列

time_threshold = 20;    % 时间间隔阈值，增大时间间隔阈值，航迹分类更准确
tolerance = 0.2;         % 容差，用于统计4秒左右的时间差
count4 = 0;              % 记录4秒左右的时间差的次数
count6 = 0;              % 记录6秒左右的时间差的次数
count8 = 0;              % 记录8秒左右的时间差的次数
count10 = 0;             % 记录10秒左右的时间差的次数
count_outrange = 0;
count_1_8 = 0;
count_dt2 = 0;
data.DeltaTime2 = NaN(height(data), 1);

% 用于保存符合条件的时间差及其所在行数
matching_deltas = {};    % 初始化为空 cell 数组

% 检查是否标签列为cell类型，如果是，则将其转换为数值
if iscell(labels)
    labels = str2double(labels);  % 转换为数值
end

labels(end+1,:) = 0;
% 如果时间是日期格式，则将其转换为秒
if isdatetime(time)
    time = seconds(time - time(1));  % 转换为秒，从第一个时间点开始
end

% 测试记录起始点和终止点索引
start_indices = find(~isnan(labels));  % 找到标签为0的索引

% 初始化存储所有时间间隔的数组
delta_ts = zeros(height(data), 1);  % 存储时间差，初始化为0

% 遍历每个测试数据
for i = 1:length(start_indices) - 1
    % 取出一个测试的所有数据 (从当前0到下一个0)
    current_start = start_indices(i);
    current_end = start_indices(i + 1) - 1;

    % 提取该段数据（按列提取）
    group_time = time(current_start:current_end);
    group_labels = labels(current_start:current_end);

    % 添加标志位，防止重复统计
    angle_exceeded = false;
    angle_exceeded_4 = false;
    angle_exceeded_6 = false;
    angle_exceeded_8 = false;
    angle_exceeded_10 = false;

    % 计算时间差
    for j = 2:length(group_time)
        delta_t = group_time(j) - group_time(j - 1);
        delta_ts(current_start + j - 1) = delta_t;  % 记录时间差

        % 检查时间间隔是否大于阈值
        if ~angle_exceeded && abs(delta_t) > time_threshold
            % 如果时间间隔大于阈值，将标签改为3
            labels(current_start + j - 1) = 3;
            count_outrange = count_outrange + 1;
            angle_exceeded = true;  % 设置标志位，避免重复统计
        end

        if ~angle_exceeded && abs(delta_t) < 1.8
            % 如果时间间隔小于1.8秒，统计
            count_1_8 = count_1_8 + 1;
            angle_exceeded = true;  % 设置标志位，避免重复统计
        end

        % 处理 DeltaTime2.时间间隔大于3秒小于20秒的轨迹有
        if ~angle_exceeded && abs(delta_t - 2) < time_threshold
            if abs(delta_t - 2) > 3
                DeltaTime2(current_start + j - 1) = 4;
                count_dt2 = count_dt2 + 1; 
                angle_exceeded = true;  % 设置标志位，避免重复统计
            end        
        end

        % 检查时间差是否在3.9到4.2秒之间，并且尚未统计过
        if ~angle_exceeded_4 && abs(delta_t - 4) < tolerance
            count4 = count4 + 1;
            angle_exceeded_4 = true;  % 设置标志位，避免重复统计
            % 保存时间差及所在行数
            matching_deltas{end + 1, 1} = current_start + j - 1;  % 行号
            matching_deltas{end, 2} = delta_t;                    % 时间差
        end

        % 类似地检查其他时间差
        if ~angle_exceeded_6 && abs(delta_t - 6) < tolerance
            count6 = count6 + 1;
            angle_exceeded_6 = true;  % 设置标志位
        end

        if ~angle_exceeded_8 && abs(delta_t - 8) < tolerance
            count8 = count8 + 1;
            angle_exceeded_8 = true;  % 设置标志位
        end

        if ~angle_exceeded_10 && abs(delta_t - 10) < tolerance
            count10 = count10 + 1;
            angle_exceeded_10 = true;  % 设置标志位
        end
    end
end

disp("测试数据超出时间间隔阈值的轨迹有：");
disp(count_outrange);
disp("测试数据的时间间隔小于1.8秒的轨迹有:");
disp(count_1_8);
disp("测试数据的时间间隔在4秒左右的轨迹有:");
disp(count4);

disp("测试数据的时间间隔在6秒左右的轨迹有:");
disp(count6);
disp("测试数据的时间间隔在8秒左右的轨迹有:");
disp(count8);
disp("测试数据的时间间隔在10秒左右的轨迹有:");
disp(count10);

disp("测试数据的时间间隔大于3秒小于20秒的轨迹有:");
disp(count_dt2);

% 将数据和标签更新回表格
labels([23170],:) = [];
labels = num2cell(labels);

data{:, 7} = labels;        % 更新标签列


% 保存修改后的表格为新的Excel文件
output_filename = 'data_test_chongfu_biao3.xlsx';
writetable(data, output_filename);

% 确认文件保存成功
disp(['数据已保存至 ', output_filename]);

然后给这些处理过的数据做标记，因为在这里我们认为含有3的数据是重新开始的一组无人机/飞鸟数据，所以我们在这里对整体的数据进行标记，其中把正常标签和标签3的数据都做标记，按顺序往下标号

%此代码用于给测试和飞鸟数据读取0到3的数据并且记录标签，序号和开始行结束行
%直接读取标签列0和3,1和3，假如不足4条则略过
clear;

% 读取Excel数据
data = readtable('data_test_chongfu_biao3.xlsx');

% 提取数据列
azimuth = data{:, 1};    % 方位角
range = data{:, 2};      % 斜距
relative_height = data{:, 3};  % 相对高度
v_h = data{:, 4};        % 径向速率
time = data{:, 5};       % 时间列
RCS = data{:, 6};        % RCS
labels = data{:, 7};     % 标签列




data.DroneNum = NaN(height(data), 1);
% 添加startline和endline列，初始为NaN
data.startline = NaN(height(data), 1);
data.endline = NaN(height(data), 1);

time_threshold = 2.5;    % 时间间隔阈值
tolerance = 0.2;         % 容差，用于统计4左右的时间差
count4 = 0;              % 记录4秒左右的时间差的次数
count6 = 0;              % 记录6秒左右的时间差的次数
count8 = 0;              % 记录8秒左右的时间差的次数
count10 = 0;             % 记录10秒左右的时间差的次数
count = 0;

droneNum = 1;
combinedData = [];
% 用于保存符合条件的时间差及其所在行数
matching_deltas = {};    % 初始化为空 cell 数组
count_31 = 0;
% 检查是否标签列为cell类型，如果是，则将其转换为数值
if iscell(labels)
    labels = str2double(labels);  % 转换为数值
end
labels(end+1,:) = 0;
% 如果时间是日期格式，则将其转换为秒
if isdatetime(time)
    time = seconds(time - time(1));  % 转换为秒，从第一个时间点开始
end

% 测试记录起始点和终止点索引
start_indices = find(~isnan(labels));  % 找到标签为0的索引
start3_indices = find(labels == 3); 
% 初始化存储所有时间间隔的数组
delta_ts = zeros(height(data), 1);  % 存储时间差，初始化为0

data.begin30 = NaN(height(data), 1);

for i = 1:length(start3_indices) - 1
     current3_start = start3_indices(i); 
     data.begin30(current3_start) = 1;
     count_31 = count_31 + 1;
end



% 遍历每个测试数据
for i = 1:length(start_indices) - 1
    % 取出一个测试的所有数据 (从当前0到下一个0)
    current_start = start_indices(i);
    current_end = start_indices(i + 1) - 1;
    
    % 找到 azimuth 列中空值的位置
    empty_rows = find(isnan(azimuth(current_start:current_end)));
    
    if ~isempty(empty_rows)
        % 如果存在空行，将 current_end 调整为空行的前一行
        current_end = current_start + empty_rows(1) - 2; % 找到第一个空行的前一行
    end
    
    azi = azimuth(current_start:current_end);


     data.DroneNum(current_start) = droneNum;
    droneNum = droneNum + 1;  % 增加无人机编号

    data.startline(current_start) = current_start;
    data.endline(current_start) = current_end;



  
       % 统计标签列的空行数
    nan_count = sum(isnan(azimuth(current_start:current_end)));

    % 标注结束行，减去空行数
    data.endline(current_start) = current_end - nan_count;

    % 提取该段数据（按列提取）
%     group_time = time(current_start:current_end);
%     group_labels = labels(current_start:current_end);
    
    % 计算时间差
%     for j = 2:length(group_time)
%         delta_t = group_time(j) - group_time(j - 1);
%         delta_ts(current_start + j - 1) = delta_t;  % 记录时间差
%         
%         % 检查时间间隔是否大于阈值
%         if abs(delta_t) > time_threshold
%             % 如果时间间隔大于阈值，将标签改为3
%             labels(current_start + j - 1) = 3;
%         end
% 
%         % 检查时间差是否在3.9到4.2秒之间
%         if abs(delta_t - 4) < tolerance
%             count4 = count4 + 1;
%             % 保存时间差及所在行数
%             matching_deltas{end + 1, 1} = current_start + j - 1;  % 行号
%             matching_deltas{end, 2} = delta_t;                    % 时间差
%         end
%         
%         % 类似地检查其他时间差
%         if abs(delta_t - 6) < tolerance
%             count6 = count6 + 1;
%         end
% 
%         if abs(delta_t - 8) < tolerance
%             count8 = count8 + 1;
%         end
%         
%         if abs(delta_t - 10) < tolerance
%             count10 = count10 + 1;
%         end
%     end
    % 在datatemp的birdNum列标注编号
%     datatemp.birdNum(current_start) = droneNum;
%     droneNum = droneNum + 1;  % 增加测试编号

    % 将 datatemp 添加到总的表格 combinedData 中
    datatemp = data(current_start:current_end, :);
    combinedData = [combinedData; datatemp];
end




disp("测试的轨迹共有:");
disp(length(start_indices))

disp("从标签3到下一个标签的轨迹有:");
disp(count_31);

disp("纯的标签轨迹占比:");
disp((length(start_indices)-count_31)/length(start_indices));


save combinedData.mat;

3.特征选取

不论是用深度学习的方法还是机器学习等方法，都需要进行特征提取，所以我们需要根据已知信息来提取特征，目前已知信息有物体的目标方位角、目标斜距、相对高度、径向速率、记录时间、RCS这几个物理量，我们首先通过观察三维空间下的无人机和飞鸟的轨迹，确定几个物理量，然后画出其二维直方图


仅仅选取这几个物理量是不够的，我们还需要选取，为了更科学的选取，我们从一些论文中选择了部分物理量作为特征值

剩下就不一一列举了。

4.模型训练

在选取好特征之后，就可以进入模型训练了，一开始我们试验了一些经典的模型例如决策树、随机森林、SVM等效果都不是特别好，后来又尝试了Adaboost

%% 训练模型（Adaboosting）
% 定义优化选项
optOptions = struct('AcquisitionFunctionName', 'expected-improvement-plus', ...
'MaxObjectiveEvaluations', 30, ... % 最大评价次数
'Kfold', 5); % 5 折交叉验证
classWeights = [majority_class_size / minority_class_size, 1]; % 权重反比于样本数量
% 构建集成分类模型，并对决策树的超参数进行优化
net = fitcensemble(p_train, t_train, ...
'Method', 'AdaBoostM1', ... % 使用 AdaBoost 方法
'Learners', 'tree', ... % 基学习器为决策树
'NumLearningCycles', 50, ... % 集成学习器的决策树数量
'OptimizeHyperparameters', 'auto', ... % 自动优化超参数
'Weights', classWeights(t_train), ... % 分配权重
'HyperparameterOptimizationOptions', optOptions); % 优化选项
% 查看模型信息
disp(net);

得到了初步的结果

后续在实际的测试中表现不佳，正确率就在70%左右，于是又去尝试了Xgboost，最后得到的效果较好，95%左右。

参考文献

[1]刘佳,徐群玉,陈唯实.无人机雷达航迹运动特征提取及组合分类方法[J].系统工程与电子技术,2023,45(10):3122-3131.
[2]王储君,万显荣,张伟民.基于不平衡数据集的无人机和飞鸟目标分类[C]//中国电子学会电波传播分会.第十七届全国电波传播年会会议论文集.武汉大学电子信息学院;,2022:4.
[3]管康萍,冯正康,马小艳,等.一种基于航迹特征的无人机与飞鸟目标雷达识别方法[J].上海航天(中英文),2024,41(01):130-136.
[4]汪浩,窦贤豪,田开严,等.基于CNN的雷达航迹分类方法[J].舰船电子对抗,2023,46(05):70-74.
[5]M. Rosamilia, A. Balleri, A. De Maio, A. Aubry and V. Carotenuto, “Radar Detection Performance Prediction Using Measured UAVs RCS Data,” in IEEE Transactions on Aerospace and Electronic Systems, vol. 59, no. 4, pp. 3550-3565, Aug. 2023.