基于三维应力对跛行奶牛蹄部参数的提取

刘彩霞，张 永，杨丽娟，黄 超，谢学虎

(内蒙古农业大学机电工程学院，内蒙古呼和浩特010018)

奶业是我国农业的重要组成部分，其发展对促进畜牧业产业升级、增加农民收入、提高国民身体素质都具有重要意义。截至2012年末，我国奶牛的存栏量已达到1 440万头，乳制品生产量达到2 545万t［1］。目前我国牛奶的人均消费水平不及发达国家的1/10，所以我国是世界上最具消费潜力的国家，奶业、奶牛养殖业必将进入新的增长阶段［2-3］。然而，奶牛疾病阻碍了奶业的健康发展，其中由肢蹄疾病所引起的牛跛行是仅次于乳房炎的第二大类疾病［4］。据统计，我国奶牛跛行病的发病率高达5．7% ～54．9%。它不仅导致牛身体状况下降，降低牛奶产量，而且跛行使奶牛空怀天数增加约3个周期。因此，早发现、早治疗是最有效的解决办法。人工监测奶牛跛行对于不断扩大规模的乳制品业越来越困难。奶牛的跛行会影响一些步态参数。跛行的奶牛相比于正常行走的奶牛行走缓慢，步幅较短且有更多的负重叠［5］。此外，跛行会减小奶牛对地面的峰值压力、平均压力、支撑时间和垂直冲量［6］，同时受影响的腿和健康的腿具有不对称性。因此，人们可以利用计算机自动提取奶牛的特征参数，判断跛行的奶牛，并对其进行跛行评分。

常用的奶牛跛行识别方法是视觉观察评分。牧场技术人员通过对牛弓背、点头以及牛步履等方面的变化观察，对牛进行运动评分。该方法耗时长，工作量大，易错过最佳的治疗时间。目前，已有一些研究尝试自动检测跛行。其中，采用传感器等采集奶牛生理及行为参数，提取与跛行有关的特征量，进而预测奶牛的早期跛行，虽没有干预牛的行为，操作简单，但数据量大，准确率较低［7-10］。牛跛行时，较为突出的特征是背部呈弓形、步态改变等。Pluk等［11-13］根据采集牛行走时的背部图像或获取蹄肢运动轨迹，判断奶牛跛行。Tasch等采用地面应力和步态参数分析法，对站立奶牛［14-16］和行走奶牛［17-26］采集地面应力数据，判断奶牛跛行。为了提高识别率，就必须解决由于时空上肢蹄序列重叠或相近而导致序列分离不准确、肢蹄间运动参数相关联等问题。国内关于小动物及人的步态分析研究较多，有关奶牛跛行自动识别方面的研究还处于起步阶段，没有较完整的资料供参考。

该研究的目的是利用薄膜三维应力测试系统自动采集奶牛自由行走时对地面的压力信息，计算压力中心(COP)，使用不同的识别准则及阈值条件，将时空重叠的支撑相精确分离。根据测取的信号，提取支撑相对应压力簇下与跛行有关的地面峰值压力、平均压力、传递能量、支撑时间、步长和对称性6个特征值，建立样本数据库，为后期步态分析进而判断跛行提供理论依据。

1 测试系统

测试系统由两块平行的三维应力板作为自动检测奶牛跛行系统的主要部分，左右两块完全对称的3D薄膜应力板加装在一块斜坡道上，压力板上铺有一层厚的橡胶垫以增加奶牛行走摩擦力。在两块应力板中间加装分隔板，高度适中，边角倒圆，防止奶牛在通过测试通道时左右肢蹄交叉串位，产生无效的肢蹄压力信号。

测试系统放置于挤奶站的出口。每头奶牛挤奶后从挤奶厅出口依次自由进入测试通道，挤奶厅出口装有射频读卡器，读入通过奶牛的ID号;测试通道两侧加装护栏，基本保证奶牛进入测试通道后沿直线行走，测试通道每次只允许一头奶牛通过。在护栏上加装光电传感器，通过加装在应力板过道前后的两个光电传感器，将奶牛从压力板开始通过到结束的一系列数据实时记录下来，同时获取奶牛的行走速度。这样的一系列数据代表一次试验。奶牛在每次试验前都要在固定时间称重。测试系统示意图见图1。

图1 测试系统结构示意图

基于采样频率测取每块压力板左右方向的受力(Fx)、前后方向的受力(Fy)和垂直方向的压力(Fz)及其对应的扭矩(Mx、My、Mz)的数据。利用MATLAB软件预处理被记录下的数据。为了剔除干扰信号，在接收计算机上编制滤波程序。由于地面3D应力信号的低频特性，滤波器设计为低通滤波。滤波器的截止频率通过分析预试验数据给出。该截止频率保证有用信号的正常通过，而高频干扰信号则被过滤掉。然后，将滤波后的支撑相进行特征提取和分析。通过完善人机界面，计算机能够实现接触面应力分布与时间曲线、支撑相压力中心(COP)的分析、峰值压力及对称性曲线显示等功能。

2 奶牛蹄肢支撑相识别与分类方法

根据牛的行走速度，每次试验每块应力板的长度允许一头牛留下0～4个支撑相。同一次试验奶牛蹄部的所有的支撑相(包括残缺的支撑相)被识别且牛蹄的序列能够被校正，最后通过计算获得支撑相的时间—空间特征值，以进一步进行步态分析，判断跛行奶牛。

2．1 压力中心的计算 当蹄部与应力板接触时，接触面覆盖的测力敏感元器件个数较多，通过数据采集可以获得多个三维应力信号。考虑到在局部坐标处力矩平衡的特点，可以由三维应力、力矩算出施力点的局部坐标值，再由局部坐标偏离全局坐标值即可给出施力点相对于全局坐标系的坐标值，构成一个压力中心。显然，一个支撑相是由若干个压力中心聚集在一起形成的压力簇，即一个完整的支撑相所对应的压力分布。压力坐标计算公式［27］如下:

式中，Zoff为压力板中心到橡胶垫表面的垂直距离;FX为沿X轴方向的力或X分力;FY为沿Y轴方向的力或Y分力;FZ为沿Z轴方向的力或Z分力。

2．2 支撑相识别及分类 创建一个蹄部滤波程序，测试系统的采样频率为440 Hz，设置训练样本的COP速度对应值是1 m/s，对8头奶牛作为试验样本进行测试。

2．2．1 独立支撑相的识别。对于奶牛单独蹄部构成的一个压力簇(支撑相)，在时间、空间上与其他压力簇没有互相干涉，称为独立支撑相。属于一个独立支撑相簇内各点的识别准则为:在一个周期内，压力中心周围点满足 FZ(0，1，2，…，n)-FZmin＞0(FZmin=10%G，G为奶牛体重，奶牛的G从挤奶厅的体重秤获得)且两个相邻的压力中心间距dmax-|dCOPn+1

dCOP(n)|＞0(dmax为一个压力簇内COP之间最大的距离，dmax根据奶牛的行走速度调整)［23］，通过这两个条件过滤，能够确保被识别的COP的轨迹最接近一个支撑相，这样一个支撑相就会被完整地分离出来，周边不符合以上条件的采样点被移除。

2．2．2 空间位置重叠而时间不重叠的支撑相识别。如果不同的奶牛在测试中蹄部落入同一个位置，那么可由奶牛的ID号识别。如果同一头奶牛在测试过程中蹄部踏入同一个位置，那么按照单独支撑相的准则识别，会出现两个支撑相合二为一的错误识别。分离空间在同一位置而时间不同的两个压力簇，识别准则需要再加一个条件—时间间隔设为Δtmin(压力簇之间最小的时间间隔，取值根据采样时间、奶牛的行走速度决定)，则落在同一个位置，这样时间差大于Δtmin的两个压力簇就被很好地分离为两个独立的支撑相。

2．2．3 空间和时间均接近(重叠)的支撑相的分离。当一头奶牛蹄部产生的支撑相序列在时间上出现重叠时，则COP重叠部分的轨迹是由第一个支撑相开始结束于最后一个支撑相。那么，蹄部序列则在时间上有重叠的压力簇，重叠压力簇的支撑见图2。在这个周期内，两只蹄同时在压力板上的垂直压力(Fzi)将大于只有一只蹄的压力，且两只蹄之间的距离D(n-1)＜Dn，D(n+1)＞Dn，表示一个向前移动的方向。这种两个支撑相在空间、时间上均出现重叠，从3个压力簇中分离的准则［23］为:

式中，参数D为压力簇间的距离。由图2可知，Fz曲线第一个支撑相的结束时间点(t1，n-1)不明，但是可以通过中间压力簇(t1n)的结束时间得到。类似的最后一个支撑相的开始时间也可以得到，即 t1，n-1=t1n，t0，n+1=t0n，然后第 n 个压力簇则被移除，同时n-1和n+1两个重叠的支撑相被分离出来［23］。

2．2．4 非支撑相分离滤除。在实际中，奶牛行走时会出现蹄部向前迈步过程中意外接触压力板的情况，导致测试系统产生无用的压力簇(支撑相)。因此，测试系统需要检测奶牛一只蹄是否踏到一个新位置。当蹄部意外接触压力板时，压力簇的位置将在奶牛下一步踏压力板时的预期范围之外。如果Dn-D(n-1)大于Da小于Db，那么第n个压力簇被舍去。在应用的参数标准值中，Da为负重叠时支撑相间的最大距离，取值39～41 cm;Db为运步时前后支撑相间(步幅)的最小距离，取值89～91 cm。

图2 前后支撑相重叠的时间—垂直压力示意

2．2．5 前后蹄部的支撑相分类。每个被识别的支撑相需要判断是属于奶牛前蹄部还是后蹄部的。如果对于一个蹄部序列，前蹄或后蹄对应的第一个支撑相或最后一个支撑相都是完整的支撑相，如果支撑相n-1和n接近，且︱Dn-D(n-1)︱小于Da，那么之前的支撑相n-1作用于前面的蹄部，第n个支撑相则作用于后面的蹄部。

2．2．6 完整支撑相识别。由于测试系统的灵敏性，奶牛蹄部从接触压力板到离开会产生大量的采样点，那么支撑相也会随之增多，即在一个周期T内会产生不完整的支撑相，且不完整支撑相的Fz小于Fzm(Fzm为识别一个完整支撑相的起始和结束时的最大压力值，Fz=10%G)。最后，进行周期调整。所有支撑相的起始时间和结束时间点都被校准，直到Fz变为0或开始增大，而不完整的支撑相只在残缺的部分持续延长。

测得的所有支撑相的压力特征值只有当COP值均在压力的传感器分布范围内是有效的。每项试验都是用来获得一个完整的蹄部序列。如果一头奶牛行走速度不在正常范围内，那么获得的试验数据不是稳态值，从而试验将被视为无利用价值而舍去。最后，利用计算机对所处理的支撑相进行周期校准调整，从而获得奶牛蹄部的最终支撑相数据，以便后续提取奶牛的运动特征值。

3 奶牛跛行的特征值的提取

三维应力测试系统通过上述程序对奶牛的支撑相识别和分类，采集奶牛各个蹄部完整的应力值、坐标值、时间等。测试系统能够自动提取与步态相关的特征参数，并加以分类，包括各个蹄部支撑相的数目、支撑时间、牛蹄垂悬时间、行走时间、速度、周期和步幅等。测试系统获得的支撑相经过滤波后，对其参数进行处理，可以自动计算出奶牛行走的一些特征值，同时系统能够以曲线和图表的方式显示出来。

3．1 确定评价奶牛跛行的支撑相特征值 奶牛行走时涉及的运动变量较多。为了更有效地分析奶牛步态，判断奶牛跛行，需要从这些运动变量中提取影响奶牛跛行的主要步态参数。健康奶牛(或人)的正常步态是指用自我感觉最自然、最舒适的姿势行走时的步态［28］。由表1可知，奶牛的蹄部运动变量(LMVs)对步态分析以及跛行判断有重要作用。该测试系统从LMVs中提取实用意义较大的垂直压力平均值(AGRF)、能量、支撑时间、步长、接触面的垂直压力峰值(PGRF)以及以上5个参数所对应的对称性。6个主要特征参数作为奶牛步态分析和跛行判断特征值。其中，当奶牛出现跛行或异常情况时，步态的整体特征结果被破坏，如AGRE的对称性被破坏(图3)。

图3 跛行奶牛与健康奶牛AGRE曲线

表1 健康奶牛7号和跛行奶牛5号各个蹄肢的运动参数

3．2 样本测试试验结果 通过对8头奶牛进行328次的试验测试验证，平均每头牛测试41次。其滤波程序获得的压力簇和支撑相数目见表2。

经兽医鉴定，其中一头ID号为5的奶牛左后蹄肢跛行评价分值达到3，而测试系统提取的特征值见表1。从表1可以看出，5号跛行奶牛左后肢的各个参数均与健康奶牛的或其健肢有明显差异。跛行奶牛的支撑相数目、PGRE、能量、AGRE相比于健康奶牛及其他3个蹄肢都小，而且其支持时间也短。通过对比健肢，可以很好地判断出跛行的奶牛，且能够确定具体跛行的蹄肢。

表2 压力簇及支撑相数量

4 结论

由于奶牛跛行对奶牛养殖场造成重要影响，步态分析的应用越来越广。结合步态分析研究理论，开发出一种对奶牛支撑相进行识别分类的方法，同时提取特征值。通过地面三维应力测试系统对行走奶牛的数据采集，结合MATLAB的滤波处理，对奶牛支撑相进行自动识别，并且提取出各个蹄部的完整支撑相，再对前后支撑相进行分类，最后对分离的支撑相参数进行提取，选取有利于判断跛行的6个特征值。该试验验证它能够很好地判断出跛行奶牛，同时确定跛腿。因此，利用测试系统对支撑相提取的特征值能够作为未来奶牛步态分析、跛行评分指标判断，从而为早期实现自动检测奶牛跛行，并且采取针对性的措施进行治疗提供良好的理论依据。

［1］北京智道顾问有限责任公司．2012-2013年中国奶牛养殖和乳制品行业分析与预测［R］．北京，2014:12-45．

［2］王俊勋．奶业发展动态与政策［J］．中国乳业，2011，117(9):2-5．

［3］李胜利．中国奶牛养殖产业发展现状及趋势［J］．中国畜牧杂志，2008(10):45-49．

［4］严作廷，王东升，张世栋，等．奶牛肢蹄病综合防治技术［J］．兽医导刊，2013(1):35-37．

［5］TELEZHENKOE，BERGSTENC．Influence of floor type on the locomotion of dairy cows［J］．Appl Anim Behav Sci，2005，93:183-197．

［6］RAJKONDAWAR PG，LIUM，DYERRM，et al．Comparison of models to identify lame cows based on gait and lesion scores，and limb movement variables［J］．JDairy Sci，2006，89:4267-4275．

［7］KAMPHUISC，FRANK E，BURKE JK，et al．Applying additive logistic regression todataderived fromsensors monitoring behavioral and physiological characteristics of dairy cows to detect lameness［J］．Journal of dairy science，2013，96(11):7043-7053．

［8］KRAMER E，CAVERO D，STAMER E，et al．Mastitis and lameness detection in dairy cows by application of fuzzy logic［J］．Livestock science，2009，125(1):92-96．

［9］ALSAAOD M，RÖMER C，KLEINMANNS J，et al．Electronic detection of lameness in dairy cows through measuring pedometric activity and lying behavior［J］．Applied animal behavior science，2012，142(3/4):134-141．

［10］DE MOL R M，ANDRÉ G，BLEUMER E JB，et al．Applicability of dayto-day variation in behavior for the automated detection of lameness in dairy cows［J］．Journal of dairy science，2013，96(6):3703-3712．

［11］PLUK A，BAHR C，POURSABERI A，et al．Automatic measurement of touch and release angles of the fetlock joint for lameness detection in dairy cattle using vision techniques［J］．Journal of dairy science，2012，95(4):1738-1748．

［12］SONG X Y，LEROY T，VRANKEN E，et al．Automatic detection of lameness in dairy cattle—Vision-based trackway analysis in cow′s locomotion［J］．Computers and electronics in agriculture，2008，64(1):39-44．

［13］POURSABERIA，BAHRC，PLUK A，et al．Real-time automatic lameness detection based on back posture extraction in dairy cattle:Shape analysis of cow with image processing techniques［J］．Computers and electronics in agriculture，2010，74(1):110-119．

［14］TASCH U，RAJKONDAWAR P G．The development of a SoftSeparatorTMfor a lameness diagnostic system［J］．Computers and electronics in agriculture，2004，44:239-245．

［15］PASTELL M，HAUTALA M，POIKALAINEN V，et al．Automatic observation of cow leg health using load sensors［J］．Comput Electron Agric，2008，62:48-53．

［16］GHOTOORLAR S M，GHAMSARI S M，NOWROUZIAN I，et al．Lameness scoring system for dairy cows using force plates and artificialintelligence［J］．Vet Rec，2012，170:126-130．

［17］LIU JB，DYER R M，NEERCHAL N K，et al．Diversity in the magnitude of hind limb unloading occurs with similar forms of lameness in dairy cows［J］．Journal of dairy research，2011，3:1-10．

［18］PASTELL M，HÄNNINEN L，DE PASSILLÉ A M，et al．Measures of weight distribution of dairy cows to detect lameness and the presence of hoof lesions［J］．JDairy Sci，2010，93(3):54-60．

［19］TELEZHENKO E，BERGSTEN C．Influence of floor type on the locomotion of dairy cows［J］．Applied animal behaviour science，2005，93(3/4):183-197．

［20］RAJKONDAWAR P G，LIU M，DYER R M，et al．Comparison of models to identify lame cows based on gait and lesion scores，and limb movement variables［J］．Journal of dairy science，2006，89(11):4267-4275．

［21］WALKERA M，PFAUT，CHANNONA，et al．Assessment of dairy cow locomotion in a commercial farm setting:The effects of walking speed on ground reaction forces and temporal and linear stride characteristics［J］．Research in veterinary science，2009，88(1):179-187．

［22］MAERTENSW，VANGEYTE J，BAERT J，et al．Development of a real time cow gait tracking and analysing tool to assess lameness using a pressure sensitive walkway:The GAITWISE system［J］．Bio Systems Engineering，2011，110(1):29-39．

［23］SKJØTH F，THORUP V M，DO NASCIMENTO O F，et al．Computerized identification and classification of stance phases as made by front or hind feet of walking cows based on 3-dimensional ground reaction forces［J］．Computers and electronics in agriculture，2013，90:7-13．

［24］VAN NUFFEL A，VANGEYTE J，MERTENS K C，et al．Exploration of measurement variation of gait variables for early lameness detection in cattle using the GAITWISE［J］．Livestock science，2013，156(1/2/3):88-95．

［25］RAJKONDAWARPG，TASCH U，LEFCOURTA M，et al．A systemfor identifying lameness in dairy cattle［J］．Appl Eng Agric，2002，18:87-96．

［26］PLUK A，OPSOMERG，VANWEYENBERGS，et al．Development of areal time cow gait tracking and analysing tool to assess lameness using a pressure sensitive walkway:the GAITWISE system［J］．Biosyst Eng，2011，110:29-39．

［27］李爱萍，胡军，刘合荣．基于三维测力台系统的人体步态分析［J］．微计算机信息，2010，26(13):226-228．

［28］胡雪艳，恽晓平，郭忠武，等．正常成人步态特征研究［J］．中国康复理论与实践，2006，10(12):855-857．