赵润茂，胡 炼,2※，罗锡文,2，唐灵茂，周 浩，杜 攀，贺 静

（1. 华南农业大学南方农业机械与装备关键技术教育部重点实验室，广州 510642；2. 南方粮油作物协同创新中心，长沙 410128）

0 引 言

农田高低不平引起农业机械行驶姿态变化，导致机具作业质量与效率下降。为实现机具位姿精准控制，水田激光平地机采用IMU（inertial measurement unit）实时检测拖拉机横滚角偏差并自动调节平地铲水平倾角，以使平地铲不受车身姿态影响而始终保持水平[1-4]；马敏等[5]采用红外光电传感器对烟草高度精确仿形，实现烟草打顶高度的精确控制以及打顶后抑芽剂的精量喷施；李君等[6]采用超声波冠形探测方法，设计了一种悬挂式电动柔性疏花机，其伺服控制系统能满足疏花机平面位置的精度要求；王松林[7]设计了检测喷杆高度的接触式传感装置，通过液压控制回路实现喷雾高度的自动调节。上述所有自动仿形系统均采用传感器测量被控对象实时状态，计算偏差、施加控制、测量反馈，形成闭环，反馈输入信号来源于传感器当前时刻测量值，系统偏差计算始终滞后于被控对象的一个运动周期[8]，影响即时控制响应。但通过对农机运动状态进行预测，根据预测信息设计控制律，可提高农具控制精度和农机作业效率[9]。

预测控制（model predictive control，MPC）的关键是预测模型，可根据对象的历史信息和未来输入，预测其未来状态或输出[10]。国内外诸多学者对此开展了大量研究。Mozaffari等[11]通过进化最小二乘算法预测混合动力轿车车速，以实现动力优化控制；Wu等[12]采用卡尔曼滤波器预测未来短时内车辆位置；Gallieri等[13]以AR（9）为预测模型设计了海藻收割机高程前馈控制系统；Koo等[14]以船体和飞机动力学方程为预测模型，实现了舰载机降落的预测控制；赵晓莉等[15]采用新陈代谢GM（1,1）模型预测玉米叶片的生长长度；张军等[16]采用灰色理论对温室系统温度的不精确性、时变性和多扰动等进行预测补偿及节能控制；此外，还有大量用于多旋翼无人机[17-19]和民用飞机[20]高程与姿态预测控制的预测模型研究。可见，只要是具有预测功能的信息集合，不论其具体表现形式如何均可作为预测模型，包括经典的状态方程和传递函数等[21]，但未见有为实现农机精准作业控制而进行的相关预测模型研究。考虑到农业机械田间运动姿态的随机、不精确、时变和多扰动等特性，精确的农机姿态模型难以获得，本文以水田平地机为研究对象，提出一种基于时间序列分析的、面向水田平地机平地铲预测控制律设计的预测模型在线辨识及其参数估计方法，并通过试验分析验证所建模型的准确性。

1 基于时间序列的水田平地机俯仰角预测模型

将平地机田间行驶时的俯仰角变化看作一个随机过程 Yt（t > 0），其在每一采样时刻的样本实现y1,y2,…,yt, …,yn构成一个俯仰角时间序列。

1.1 俯仰角时序模型结构识别

为识别描述平地机田间行驶俯仰角的时序模型是自回归（auto-regressive，AR）、滑动平均（moving-average，MA）或自回归滑动平均（auto-regressive and movingaverage，ARMA）模型，分别计算俯仰角数据样本的自相关函数（autocorrelation function，ACF）和偏自相关函数（partial autocorrelation function，PACF），根据其呈现的“拖尾”或“截尾”特征进行判别[22]，如表1所示。

表1 时间序列ACF和PACF的一般特征Table 1 ACF and PACF characteristics for time series

ACF和PACF的计算式[23]分别为式（1）和式（2）。式中n为样本数据量，y为样本均值，k为迟滞阶数，ρˆk为俯仰角序列ACF的估计，φˆkk为俯仰角序列PACF的估计，

由表 1和田间预试验可知，平地机田间行驶俯仰角时序的 ACF和 PACF均呈现“拖尾”现象，因此采用ARMA表达平地机俯仰角动态变化，模型结构为式（3）。

式中 A （ q ） = 1 + a1q-1+…+，

C（ q ） = 1 + c1q-1+…+q-nc，y（ t）为t时刻模型输出，na为自回归项阶数，nc为滑动平均项系数，q为后移算子， e（ t）为系统噪声， e（ t）~N ID(0， σ2)。

1.2 俯仰角ARMA模型定阶

采用基于信息量准则的 AIC（akaike information criterion）定阶方法[24]估计平地机俯仰角预测模型ARMA的阶次（na, nc），式（4）。

式中N为样本个数；为残差方差。

AIC定阶需要进行大范围内的阶数搜索，并对每次搜索得到的模型计算AIC信息量[25]；随着搜索阶数的递增，式（4）中第一项数值下降，第二项数值上升，直至某阶AIC（na, nc）最小即为最佳阶数。

1.3 模型参数在线辨识

模型结构确定后，需根据平地机俯仰角传感信息估计式（3）中的参数 a1, a2, …,c1, c2,… cnc值。水田地况复杂多变，采用批处理最小二乘法（least squares，LS）离线辨识得到的定常参数模型无法长期描述平地机俯仰角的动态变化；而始终利用一定长度的最新数据进行在线参数估计虽可有效跟踪系统时变，但运算速度无法满足预测控制实时性要求，尤其是ARMA中MA部分的参数估计需要大量非线性计算[26]。为保证平地机俯仰角的时变模型参数实时有效估计又避免每次传感数据更新均需重复估计参数的耗时运算[27]，采用平地机最新俯仰角数据对上一步估计出的参数进行校正，即递推最小二乘算法（recursive least square，RLS）[28]。RLS参数递推方程为式（5）～（7）。

式中( t) = [ a1, … ,,c1, … cn]为t时刻模型拟辨识参数；c

yˆ( t)为t时刻平地机俯仰角估计值； y( t)为t时刻平地机俯仰角测量值；φˆ（ t ） = [- y （ t- 1 ）,… ,- y （ t- na）,（ t - 1）,… ,eˆ（ t - nc）]T， eˆ ( t) = y( t) -( tˆ(t )，为t时刻模型残差；K( t)为校正系数。

2 水田平地机俯仰角预测模型构建

2.1 数据采集与预处理

试验在华南农业大学岑村实验基地进行，将 AHRS（XSENS, Mti-300）安装于以久保田2ZGQ-6D1为配套动力的1PJ-3.0型水田激光平地机座椅后方，采集平地机行驶俯仰角数据，采样频率为20 Hz。传感数据通过串口通信RS232发至上位机实时处理，平地机俯仰角在线建模在MATLAB 2016b平台中实现。

初始采集得到的数据由于可能包含趋势性或周期性等因素，无法保证时间序列的统计特征不随时间推移而变化[29]，故不能直接用于时序建模。因此，采用差分法去除平地机俯仰角时序中非平稳因素，建模后通过逆差分还原序列[30]。图1为任意采集的一段时长为200 s的平地机水田行驶俯仰角度值，由图可见，在水田中平地机行驶姿态起伏变化频繁，为非平稳随机过程，通过一阶差分使序列平稳化，差分处理后俯仰角数据围绕零均值上下波动，序列基本平稳。序列平稳性采用MATLAB提供的检验函数adftest进行ADF单位根检验，函数返回值为1表明序列平稳，可以进行时序建模。

图1 平地机水田行驶俯仰角序列Fig.1 Time series for pitch angle of leveler in paddy field

2.2 预测模型结构识别

差分后序列ACF与PACF值的计算结果如图2所示。图2表明，俯仰角序列的ACF与PACF值绝大部分落在2倍标准差外，其衰减过程震荡且收敛缓慢，可确定ACF与 PACF“拖尾”，依据表 1，选择 ARMA（na，nc）（na>0，nc>0）模型对平地机俯仰角序列建模[23]。从na= 1，nc= 1开始逐次升阶拟合前20阶平地机俯仰角的ARMA（na，nc）模型并计算AIC值。对于高阶模型定阶，选用（2n，2n–1）的建模方案[26]进行最佳阶数搜索以减小运算次数，AIC计算结果如表2所示，最小AIC值为–5.272 9，此时模型阶数na=18，nc=17。再在（18,17）附近各进行一步搜索（19,18）与（17,16）模型，计算结果AIC（19，18）=–5.269 5，AIC（17，16）= –5.271 9。综上可确定最佳模型为：ARMA（18，17）。

图2 差分序列ACF与PACFFig.2 ACF and PACF of differential time series

2.3 预测模型参数估计

采用RLS算法在线辨识模型参数，令式（5）～（7）中初值ˆ（0）=0，P(0)=μI，其中I为35阶单位矩阵，μ=1× 1 06；可以证明当μ取足够大的正数时，经过不多次的递推，初值影响就会消失[27]。输出俯仰角 ARMA（18,17）模型的向前1步递推值并与AHRS实测值比较，如图3所示，误差分析表明最大绝对误差MAE（maximum absolute error）为 3.527 1°，均方根误差 RMSE（root mean square error）为 1.274 9°。

表2 AIC计算结果Table 2 Results of AIC calculation

图3 RLS辨识模型输出结果Fig.3 Output results of model identified by RLS

对图3的分析可看出，约50 s后模型绝对误差开始增大，在近150 s处达到最大为–3.527 1°，这是因为随着数据量的增加，P （ t）、K （ t）变小，产生“数据饱和”现象，从而对θˆ（t）的修正能力变弱，导致时变参数估计失败。为改善RLS在时变系统中的跟踪性能，提出采用遗忘因子[28,31]（FF，forgetting factor）λ（0＜λ≤1）对在线传感数据加权以强调新息对参数估计的贡献，保证RLS可跟踪系统参数变化。FFRLS辨识算法见式（8）～式（10）。

FFRLS算法中遗忘因子λ的取值对俯仰角模型精度有重要影响。若λ选取过小，模型过多重视新息，模型输出较实测值跳动震荡大；若λ选取过大，即建模更注重序列历史长期趋势，模型输出则较为平滑，对平地机俯仰角值的高频变化表达不够。一般地，若时间序列的连续T0个采样趋于稳定常值，则有式（11）[32]。

对于时变与非线性较强的平地机俯仰角序列，T0=3，此时0.67λ=。为更加精确确定λ值，分别在0.61λ≤≤时进行仿真试验并以模型误差RMSE与MAE评价λ取值优劣，结果如表3所示。依据表3可见，对于当前序列，当遗忘因子0.7λ=时，模型输出RMSE为0.026 1°，MAE为0.027 0°，误差最小。

2.4 预测模型适用性检验

在将模型用于平地机俯仰角预测前，应对其适用性进行检验，主要包括残差分析与模型仿真。残差分析检验残差是否为白噪声，若为白噪声，说明姿态模型拟合良好，否则说明残差序列仍有相关信息未被提取，需要重新建立模型[33]；通过模型仿真可对模型输出误差进行直观分析与比较。根据式（1），ARMA模型的1步递推输出为式（12）所示。

表3 λ取值对模型拟合残差的影响Table 3 Influence of different λ values on residual error

写成最小二乘形式为式（13）所示。

对ARMA（18，17）模型进行残差自相关检验，如图4所示。图4表明高于1阶滞后的残差自相关函数值均落在非显著相关的99%置信区间内，说明残差序列值之间不存在任何相关性，为白噪声，因此当前 ARMA（18，17）模型适用于平地机俯仰角序列建模。

图4 残差相关性检验Fig.4 Correlation test of residual

3 试验与结果分析

为进一步验证 ARMA（18，17）模型对描述平地机俯仰角的普适性与鲁棒性，于2017年5月31日在华南农业大学岑村实验基地（北纬 23.163429°，东经113.370198°）任选取包含坑洼、坡度、路障的自然水泥路面以模拟阶跃、斜坡、正弦激励信号源。驾驶 1PJ-3.0型水田激光平地以作业速度（1.2 m/s）分别完成过沟，爬坡，前进-后退越障、前进-掉头-前进越障，每组试验重复3次，如图5所示。利用AHRS实时采集上述地况下平地机行驶俯仰角度值并同步在线建模，采样率为 20 Hz，依次比较模型输出与AHRS实测值，如图6a～6d所示，图6e为平地机水田（旋耕后静置48 h）作业时连续200 s俯仰角测量值与模型输出对比。

从图6看出，在4种模拟典型地况及水田环境下，算法的一步递推输出结果与AHRS实测值趋势变化相同，所建模型可有效跟踪平地机俯仰角动态变化，具有较高准确性。误差分析结果如表 4所示，各模拟地况下（图6a～图6d）模型输出MAE与RMSE值均小于0.1°；水田作业时(图6e)模型输出MAE与RMSE值小于0.2°，误差较大是由于复杂环境中模型时变性强造成。综上所述，用ARMA（18，17）对平地机俯仰角建模并通过FFRLS实时估计和递推更新模型参数的方法是可行的，所建模型输出精度满足平地机姿态预测控制律设计的要求。

图5 不同地况模拟试验Fig.5 Algorithm testing under different situations

图6 不同地况下ARMA（18，17）预测与AHRS实测对比Fig.6 Comparison between output value of ARMA(18, 17)and AHRS measurement

表4 不同模拟地况下ARMA（18,17）输出误差Table 4 Output errors of ARMA (18,17) under different situations

4 结 论

1）基于时间序列分析方法，建立了水田平地机俯仰角序列的 ARMA（18，17）模型，通过模型拟合残差单位根ADF检验，验证了所建模型适用性。

2）为满足预测控制系统对预测模型辨识的实时性要求，采用遗忘因子递推最小二乘法对模型参数在线估计并递推更新，通过模型仿真试验得到其遗忘因子取值0.7。分别采集几种模拟典型地况及实际水田作业时平地机行驶俯仰角数据在线建模，将模型输出值与 AHRS实测值比较，结果表明：模型输出与AHRS实测值变化趋势一致，MAE与RMSE计算值均不超过0.2°，说明采用ARMA（18，17）对平地机俯仰角建模并通过FFRLS在线估计模型参数的方法是可行的，满足控制系统对平地机俯仰角预测的误差要求。

面向水田平地机姿态预测控制律设计，本文提出了平地机俯仰角预测模型在线辨识方法。由于水田平地过程中，平地铲运动并不影响车身俯仰角变化，因此本文研究平地机俯仰角在线建模时，未分析平地铲与动力车身之间的相对运动关系。后续在设计MPC控制器时，将针对性开展铲姿态与车身俯仰角之间的运动关联分析，并结合硬件平台测试系统动态响应，进一步对建模算法精度和执行效率进行检验、分析与优化。

[1] 胡炼，林潮兴，罗锡文，等. 农机具自动调平控制系统设计与试验[J]. 农业工程学报，2015，31(8)：15－20.Hu Lian, Lin Chaoxing, Luo Xiwen, et al. Design and experiment on auto leveling control system of agricultural implements[J].Transactions of the Chinese Society of Agricultural Engineering(Transactions of the CSAE), 2015, 31(8): 15－20. (in Chinese with English abstract)

[2] 胡炼，罗锡文，林潮兴，等. 1PJ-4.0型水田激光平地机设计与试验[J]. 农业机械学报，2014，45(4)：146－151.Hu Lian, Luo Xiwen, Lin Chaoxing, et al. Development of 1PJ-4.0 Laser Leveler Installed on a Wheeled Tractor for Paddy Field[J]. Transactions of the Chinese Society for Agricultural Machinery, 2014, 45(4): 146－151. (in Chinese with English abstract)

[3] 赵祚喜，罗锡文，李庆，等. 基于MEMS惯性传感器融合的水田激光平地机水平控制系统[J]. 农业工程学报，2008，24(6)：119－124.Zhao Zuoxi, Luo Xiwen, Li Qing, et al. Leveling control system of laser-controlled land leveler for paddy field based on MEMS inertial sensor fusion[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE),2008, 24(6): 119－124. (in Chinese with English abstract)

[4] 胡炼，罗锡文，林潮兴，等. 一种作业机具自动调平系统及其控制方法：CN104206063A[P]. 2014-12-17.

[5] 马敏，张晓辉，宋涛，等. 智能烟草打顶抑芽机控制系统[J]. 农业工程学报，2011，27(2)：129－135.Ma Min, Zhang Xiaohui, Song Tao, et al. Control system of intelligentized tobacco chopping and sucker controlling machine[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2011, 27(2): 129－135. (in Chinese with English abstract)

[6] 李君，徐岩，许绩彤，等. 悬挂式电动柔性疏花机控制系统设计与试验[J]. 农业工程学报，2016，32(18)：61－66.Li Jun, Xu Yan, Xu Jitong, et al. Design and experiment of control system for suspended electric flexible thinner[J].Transactions of the Chinese Society of Agricultural Engineering(Transactions of the CSAE), 2016, 32(18): 61－66. (in Chinese with English abstract)

[7] 王松林. 喷杆高度自动调节系统设计与试验研究[D]. 杨凌：西北农林科技大学，2014.Wang Songlin. Design and Experiments on Boom Height Automatic Adjusting System[D]. Yangling: Northwest A&F University, 2014. (in Chinese with English abstract)

[8] 林守利. 随机海浪作用下船上物体保持水平的预测控制研究[D]. 沈阳：沈阳工业大学，2003.Lin Shouli. The Study on Predictive Control of an Object's Level on a Ship in Random Ocean Waves[D]. Shenyang:Shenyang University of Technology, 2003. (in Chinese with English abstract)

[9] Suh J, Yi K, Jung J, et al. Design and evaluation of a model predictive vehicle control algorithm for automated driving using a vehicle traffic simulator[J]. Control Engineering Practice, 2016, 51: 92－107.

[10] 席裕庚. 预测控制(第2版)[M]. 北京：国防工业出版社，2013.

[11] Mozaffari L, Mozaffari A, Azad N L. Vehicle speed prediction via a sliding-window time series analysis and an evolutionary least learning machine: A case study on San Francisco urban roads[J]. Engineering Science and Technology, an International Journal, 2015, 18(2): 150－162.

[12] Wu Chaozhong, Peng Liqun, Huang Zhen, et al. A method of vehicle motion prediction and collision risk assessment with a simulated vehicular cyber physical system[J]. Transportation Research Part C: Emerging Technologies, 2014, 47: 179－191.

[13] Gallieri M, Ringwood J, Giantomassi A, et al. Predictor design for altitude control of a seaweed harvester[J]. IFAC Conference on Control Applications in Marine Systems,Germany, 2010, 43(20): 301－306.

[14] Koo S, Kim S, Suk J. Model predictive control for UAV automatic landing on moving carrier deck with heave motion[J]. IFAC Proceedings Volumes (IFAC-Papers Online), 48(5): 59－64.

[15] 赵晓莉，马新明，王举才，等. 基于新陈代谢 GM(1，1)模型的玉米叶长动态预测[J]. 农业工程学报，2013，29(10)：183－189.Zhao Xiaoli, Ma Xinming, Wang Jucai, et al. Dynamic prediction on leaf length of maize based on metabolic model GM(1, 1)[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2013, 29(10):183－189. (in Chinese with English abstract)

[16] 张军，张侃谕. 温室温度控制系统不确定性与干扰的灰色预测补偿算法[J]. 农业工程学报，2013，29(10)：225－233.Zhang Jun, Zhang Kanyu. Grey prediction compensation algorithm for the uncertainty and interference of greenhouse temperature control system[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE),2013, 29(10): 225－233. (in Chinese with English abstract)

[17] Bangura M, Mahony R. Real-time model predictive control for quadrotors[J]. IFAC Proceedings Volumes, 2014: 11773－11780.

[18] Alexis K, Nikolakopoulos G, Koveos Y, et al. Switching model predictive control for a quadrotor helicopter under severe environmental flight conditions[J]. IFAC Proceedings Volumes, 2011, 44(1): 11913－11918.

[19] Alexis K, Nikolakopoulos G, Tzes A. Switching model predictive attitude control for a quadrotor helicopter subject to atmospheric disturbances[J]. Control Engineering Practice,2011, 19: 1195－1207.

[20] Dong Yiqun, Zhang Youmin, Ai Jianliang. Full-altitude attitude angles envelope and model predictive control-based attitude angles protection for civil aircraft[J]. Aerospace Science and Technology, 2016, 55: 292－306.

[21] 邹涛，丁宝苍，张端. 模型预测控制工程应用导论[M]. 北京：化学工业出版社，2010.

[22] Cryer, Jonathan D, Chan, et al. Time Series Analysis: with Applications in R[M]. Berlin: Springer Science & Business Media, 2008.

[23] 王燕. 应用时间序列分析[M]. 北京：中国人民大学出版社，2015.

[24] 彭秀艳，赵希人，高奇峰. 船舶姿态运动实时预报算法研究[J]. 系统仿真学报，2007，19(2)：267－271.Peng Xiuyan, Zhao Xiren, Gao Qifeng. Research on real-time prediction algorithm of ship attitude motion[J]. Journal of System Simulation, 2007, 19(2): 267－271. (in Chinese with English abstract)

[25] 龚啸. 基于 ARMA递推算法的电力系统低频振荡模式在线辨识研究[D]. 重庆：重庆大学，2011.Gong Xiao. Research on Online Identification for Low Frequency Oscillation Modes in Power System Based on ARMA Recursive Method[D]. Chongqing: Chongqing University, 2011. (in Chinese with English abstract)

[26] 潘迪特，吴宪民. 时间序列及系统分析与应用[M]. 北京：机械工业出版社，1988.

[27] 高奇峰. 舰载机起降指导系统仿真研究[D]. 哈尔滨：哈尔滨工程大学，2006.Gao Qifeng. The Simulated Research for Guidance System of Carrier Plane Taking Off & Landing[D]. Harbin: Harbin Engineering University, 2006. (in Chinese with English abstract)

[28] Ljung L. System identification：theory for the user[M].Beijing: Tisinghua University Press, 2002.

[29] 钟美霞. 不同采样频率时间序列的多步预测比较研究[D].广州：广东工业大学，2013.Zhong Meixia. A Comparative Study of the Multi-step Prediction of Time Series with Different Sampling Frequencies[D].Guangzhou: Guangdong University of Technology, 2013. (in Chinese with English abstract)

[30] 杨叔子，吴雅，轩建平. 时间序列分析的工程应用[M]. 武汉：华中科技大学出版社，2007.

[31] 庞中华，崔红. 系统辨识与自适应控制[M]. 北京：北京航空航天大学出版社，2013.

[32] Mathworks. Matlab Product Documentation[Z]. Natick: 2016.

[33] 苏孟辉. 基于时间序列的网络流量监测系统的设计与实现[D]. 广州：华南理工大学， 2015.Su Menghui. Design and Implementation of Network Traffic Monitoring System Based on Time Series[D]. Guangzhou:South China University of Technology, 2015. (in Chinese with English abstract)