CN109540089B - 一种基于贝叶斯-克里金模型的桥面高程拟合方法 - Google Patents

一种基于贝叶斯-克里金模型的桥面高程拟合方法 Download PDF

Info

Publication number
CN109540089B
CN109540089B CN201811200557.4A CN201811200557A CN109540089B CN 109540089 B CN109540089 B CN 109540089B CN 201811200557 A CN201811200557 A CN 201811200557A CN 109540089 B CN109540089 B CN 109540089B
Authority
CN
China
Prior art keywords
model
bayesian
kriging
fitting
elevation
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201811200557.4A
Other languages
English (en)
Other versions
CN109540089A (zh
Inventor
余晓琳
贾布裕
颜全胜
杨铮
陈宇轩
杨钰炜
罗宇蕃
黄逸锋
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
South China University of Technology SCUT
Original Assignee
South China University of Technology SCUT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by South China University of Technology SCUT filed Critical South China University of Technology SCUT
Priority to CN201811200557.4A priority Critical patent/CN109540089B/zh
Publication of CN109540089A publication Critical patent/CN109540089A/zh
Application granted granted Critical
Publication of CN109540089B publication Critical patent/CN109540089B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C5/00Measuring height; Measuring distances transverse to line of sight; Levelling between separated points; Surveyors' levels
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/10Geometric CAD
    • G06F30/13Architectural design, e.g. computer-aided architectural design [CAAD] related to design of buildings, bridges, landscapes, production plants or roads
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06NCOMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00Computing arrangements based on biological models
    • G06N3/02Neural networks
    • G06N3/08Learning methods
    • G06N3/084Backpropagation, e.g. using gradient descent

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Geometry (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • General Engineering & Computer Science (AREA)
  • Computational Linguistics (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Molecular Biology (AREA)
  • Computing Systems (AREA)
  • Biophysics (AREA)
  • Biomedical Technology (AREA)
  • Mathematical Physics (AREA)
  • Software Systems (AREA)
  • Artificial Intelligence (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Data Mining & Analysis (AREA)
  • Architecture (AREA)
  • Civil Engineering (AREA)
  • Structural Engineering (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Complex Calculations (AREA)

Abstract

本发明桥梁工程中桥面高程测量技术领域,公开了一种基于贝叶斯‑克里金模型的桥面高程拟合方法,包括以下步骤:S1、建立贝叶斯‑克里金拟合模型;S2、高程测试点样本试验优化;S3、建立贝叶斯‑克里金预测模型;S4、对全桥面进行高程拟合评估。其有益效果在于,将贝叶斯和克里金结合:基于非信息超先验,对克里金模型的基函数系数、相关参数赋予了多层先验约束,利用EM算法,求解基函数系数、相关参数的最大后验估计,对克里金模型进行改进,建立贝叶斯‑克里金模型,增强了模型的自适应性和稳健性。同时进行基于数论法的高程测量样本的试验优化设计。

Description

一种基于贝叶斯-克里金模型的桥面高程拟合方法
技术领域
本发明涉及桥梁工程中桥面高程测量技术领域,具体涉及一种基于贝叶斯-克里金模型的桥面高程拟合方法。
背景技术
桥梁在使用过程中,在内部因素(如收缩徐变,材料老化等)以及外部因素(如交通荷载,温度作用等)作用下,桥面高程会发生变化(沉降或上挠)。对桥梁桥面挠度或高程进行检测时,不可能做到对全桥面每个点进行测量,只能选取小数局部点(如跨中、4分跨、8分跨等)进行测量。事实上由于桥梁在内外部因素的复杂影响作用下,其桥面高程分布是不均匀分布的,那些局部高程点数据信息并不能有效描述全桥面的实际状况,尤其在评估其修补需要的方量(很可能得到与实际相差较大的结果)。
现有技术中一般采用以下几种桥面高程拟合的方法:
1、一元线性回归分析法:
一元回归模型是针对单一自变量和单一因变量间的近似线性关系,通过一元方程完成拟合,进而运用获得的一元线性方程去获得预计结果。在基于一元回归分析的桥面高程拟合方法中,一般选择桥梁纵向坐标作为自变量x,竖向坐标作为因变量y,不考虑横向分布,这是一种最简单的处理方式,其原理如下:
建立一元线性回归的数学模型:
y=a+bx+Δ (1-1)
式中,y-竖向坐标(高程);x-桥梁纵向坐标;a,b-回归系数;Δ-随机误差。
已知n对局部测试点数据:{(x1,y1);(x2,y2);...(xn,yn)},可采用最小二乘法对回归系数a,b进行求解估计,从中可以得到预测模型的拟合曲线以及系数估计值(最大似然估计值)。
2、二元线性回归分析法:
在基于二元回归分析的桥面高程拟合方法中,考虑横向分布,选择桥梁纵向坐标作为自变量x,桥梁横向坐标作为自变量y,竖向坐标作为因变量z。
其模型为
z=a+bx+cy+Δ) (1-2)
根据n对局部测试点数据:{(x1,y1,z1);(x2,y2,z2);...(xn,yn,zn)},并结合最小二乘法,可求得系数a,b,c,进而建立预测模型。
3、趋势面分析法:
把实际的地理曲面分解为趋势面和剩余面两部分,前者反应空间位置的宏观分布规律,属于确定性因素作用的结果;而后者则对应于微观区域,被认为是随机因素影响的结果。趋势面分析的一个基本要求就是,所选择的趋势面模型应该是剩余值最小,而趋势值最大,这样拟合度精确度才能达到足够的准确性。趋势面分析是通过回归分析原理,运用最小二乘法拟合一个二维非线性函数,模拟空间要素在空间上的分布规律。趋势面分析的核心就是从实际观测值出发推算趋势面,一般采用回归分析方法,使得残差平方和最小从而估计趋势面参数。这里以二次趋势面为例,简单介绍其原理:
设高程z和桥梁纵向坐标x以及桥梁横向坐标y,有以下关系:
z=a0+a1x+a2y+a3x2+a4y2+a5xy+Δ (1-3)
根据n对局部测试点数据:{(x1,y1,z1);(x2,y2,z2);...(xn,yn,zn)},将预测模型表达式写成矩阵形式,
Z=XA+Δ
Figure BDA0001829890230000021
因为未知系数有6个,因此至少需要6个已知测量点数据才能求解此方程。根据最小二乘原理可求待定参数:A=(XTX)-1XTZ。建立好预测模型后,可利用该模型预测桥面其他点位的高程。
4、最近邻点法:
又称之为泰森多边形方法。它采用一种极端的边界内插方法-只用最近的单个点进行区域插值(区域赋值)。该方法按数据点位置将区域分割成子区域,每个子区域包含一个数据点,各子区域到其内数据点的距离小于任何到其它数据点的距离,并用其内数据点进行赋值。其模型公式主要可写为:
ze=zi (1-5)
ze为待估点值,zi为i测点值。i点须满足:dei=min(de1,de2,...,den),其中dei为空间两点的欧几里得距离。
5、样条函数插值法:
利用一种灵活的分段曲线逐段的拟合出平滑的曲线。这种分段曲线称为样条。曲线规绘出的曲线在数学上用分段的三次多项式函数来描述这种曲线,其连接处有连续的一阶和二阶连续导数。一般的分段多项式p(x)定义为:
Figure BDA0001829890230000031
x1,...,xk-1将区间x0,xk分成k个子区间,这些分割点称“断点”,曲线上具有这些x值的点称为“节”。函数pi(x)为小于等于m次的多项式。r项用来表示样条函数的约束条件:r=0,无约束;r=1,函数连续且对它的导数无任何约束;r=m-1,区间[x0 xk]可用一个多项式表示;r=m,约束条件最多。当r=m=3,分段曲线是三次多项式,因此又被人们称为三次样条函数。在桥面高程拟合中,可采用三次样条函数进行拟合。
6、BP神经网络法:
BP(Back Propagation)神经网络网络是一种按误差逆传播算法训练的多层前馈网络。BP神经网络结构通常分为输入层、隐含层和输出层,学习过程主要分为正向和反向传播两个阶段:信号的正向传播阶段和信号的反向传播阶段。在桥面高程拟合计算,其主要流程如下:
第一步:网络初始化,初始化权值和闽值。设定阀值的初值,并即对所有权值Wij赋以小的随机数;设定纵向坐标x和横向坐标y为输入量,竖向坐标z(高程)为输出量;
第二步:给定训练数据集{(x1,y1,z1);(x2,y2,z2);...(xn,yn,zn)},提供输入向量qi=(xi,yi)及其对应的期望输出zi
第三步:计算实际输出zi=f(∑Wijqi)。式中,f为标准Sigmoi函数,i和j分别表示输人节点和隐含节点;
第四步:调整权值,按照误差反向传播方向,从输出节点开始将误差反馈到隐含层,并按照下式对权值Wij修改:Wij(t+1)=Wij(t)+γδjzj
式中:γ为大于0的增益;δj为对应节点J的误差;t为迭代次数。
第五步:判断网络误差是否达到预定精度。当网络的实际输出达到预定的精度或者达到预先设定的学习次数,算法结束。否则,返回到第二步,进入下一轮学习。
在桥面高程拟合除了上述几种主要拟合方法外,还可以采用移动平均插值方法,径向基函数插值法,最小曲率法等。
以上几种现有技术中存在的问题:
1、一元线性回归分析法、二线性回归分析法以及趋势面分析法实质上属于整体插值法,需要使用方差分析和回归方程等标准的统计方法,对数据要求高,计算量较大,结果点很少通过原始数据点,只是对整个研究区产生最佳拟合效果。其中,一元线性回归分析法不能考虑横向坐标的影响,其拟合结果最为粗糙,但其原理方法简单(可直接利用EXCEL处理),在实际工程中经常被使用。趋势面分析法效果取决于实际问题,高次趋势面的效果并不一定比低次好,而且高趋势面由于含有高次项,容易在在数据区外围产生异常高值或低值。除此之外,这些模型还存在一个重要问题就是这些方法的结构形式较为刻板不具备自适应性,对于不同的问题应该选取什么样的参数表达式,不同的参数表达式对拟合的计算结果有多大影响是不清楚的。
2、最近邻点法、样条函数插值法实质上属于局部插值拟合技术,其中,最近邻点法受样本点(已有测试点)的影响较大,只考虑距离因素,对其他空间因素和变量所固有的某些规律没有过多地考虑。实际应用中,效果常不十分理想;而样条函数插值法不适用于在短距离内属性有较大变化的地区,否则估计结果偏大;不适用于在短距离内属性有较大变化的地区,否则估计结果偏大。他们同样会存在和上述方法一样的问题:形式较为刻板不具备自适应性,无法有效知道参数对拟合模型的影响作用。
3、BP神经网络法属于机器学习范畴,但神经网络并没有严格的数学基础,存在以下问题:需要预先设定神经网络的结构或在训练过程中不断地进行摸索,而使这种方法对“使用者”先验知识和经验过分的依赖,造成结构选择困难问题;BP神经网络可能会陷入局部极小,造成局部极值问题;BP神经网络对训练数据的需求很大,需要较多测试数据的支持,但在实际桥面拟合测试中,测试获得的数据是较为有限的,这时很可能会产生“欠学习”问题。
发明内容
本发明的目的是为了克服以上现有技术存在的不足,提供了一种利用较少数量的桥面高程测点,可对全桥范围内的其他桥面高程点进行高精度拟合的基于贝叶斯-克里金模型的桥面高程拟合方法。
本发明的目的通过以下的技术方案实现:一种基于贝叶斯-克里金模型的桥面高程拟合方法,包括以下步骤:
S1、建立贝叶斯-克里金拟合模型;
S2、高程测试点样本试验优化;
S3、建立贝叶斯-克里金预测模型;
S4、对全桥面进行高程拟合评估。
进一步地,所述S1包括以下内容,
S11、建立贝叶斯-克里金模型回归系数求解子系统:
S111、以范-克里金模型y(x)=f(x)Tβ+z(x)为基础,式中:f(x)=[f1(x)…fp(x)]表示已经选择的基函数,β为回归系数,z(x)是服从均质为0、协方差矩阵为σ2R的高斯过程;得到y*=F*β+z*,此时z*是均值为0、协方差矩阵为σ2I单位矩阵的高斯过程;
S112、基于贝叶斯三层先验的克里金回归系数建模:
第一层为对随机过程的方差逆(1/σ2)赋予Gamma先验:p(a|u,v)=Γ(a|u,v),式中:a=1/σ2,u,v是Gamma先验中的超参数;
对回归系数β赋予变方差的高斯先验:
Figure BDA0001829890230000053
式中:N(βi|0,τi)是每一个βi对应的一个满足均值为零,参数τi为高斯密度函数,
Figure BDA0001829890230000051
第二层为对参数τi赋予Gamma先验:
p(τi|λ)=Γ(τi|1,λ/2)=0.5λexp(-0.5λτi),
第三层为对参数λ赋予Gamma先验:
p(λ|η)=Γ(λ|0.5η,0.5η),式中:η表示超参数;
得出β的估值为:
Figure BDA0001829890230000052
S12、建立贝叶斯-克里金模型相关参数求解子系统:
S121、以范-克里金模型y(x)=f(x)Tβ+z(x)为基础,输出y的似然函数为:
Figure BDA0001829890230000061
S122、基于贝叶斯两层先验的克里金相关参数建模:
第一层为相关参数θ赋予变方差的高斯分布;
首先对回归系数β和方差逆a赋予如下常数先验:
Figure BDA0001829890230000062
为了得到自适应解,对相关参数θ提出一个多层先验:
Figure BDA0001829890230000063
式中:N(θi|0,ωi)是一个变方差的高斯分布,ωi表示方差,
Figure BDA0001829890230000064
第二层为对每一个ωi取Jeffreys非信息的超先验:
p(ωi)=1/ωi
进一步地,所述S2中采用数论方法对高程测试点样本试验进行优化。
进一步地,所述S2中采用正交试验、均匀设计、中心复合设计或拉丁超立方抽样中的任一种方法对高程测试点样本试验进行优化。
进一步地,所述S4包括以下内容:
向S3的贝叶斯-克里金预测模型中,输入其他桥面点的平面位置,拟合估算出对应的高程。
进一步地,所述S4还包括利用S3的贝叶斯-克里金预测模型预测均方差对拟合效果进行评估。
本发明的发明原理为:
1)基于贝叶斯-克里金模型的桥面高程拟合技术
已有技术在对桥面高程进行拟合时,以多项式回归法和样条插值法为主,其中多项式回归法利用已有的高程点数据进行最小二乘求解,得到一个代表全局分布的近似拟合表达式,其拟合精度取决于测点数据,并不注重局部区域的拟合。而样条插值法采用分段曲线逐段拟合出平滑的曲线,其注重考虑距离因素,对其他空间因素和变量所固有的某些规律没有过多地考虑,即强调局部区域的拟合,在短距离内属性有较小变化的地区拟合效果较好。无论是多项式回归法还是样条插值法,他们的拟合效果对数据样本的选取非常敏感,当存在异常点时(如个别点高程过低或过高),他们的拟合效果都将受到影响;多项式回归法对数据要求较为苛刻,它的拟合精度取决于数据样本的质量而非数量,而样条插值法却需要大数据的支持,一般情况下数据信息越大,其拟合效果越好。
本发明技术对桥面高程进行拟合时,首先利用了克里金模型优秀的预测值估计性能及预测均方差估计(可用于评估预测效果好坏)能力,克里金模型是一种基于随机过程的统计预测法,具有平滑效应以及估计方差最小的统计特征,被视为最优的线性无偏估计。克里金模型作为线性回归分析的一种改进的技术,包含了线性回归部分和非参数部分,其中线性回归部分提供模拟的全局近似,而非参数部分被视作随机分布的实现,提供模拟的局部近似,是一种全局和局部相对平衡的拟合技术,在小样本的情况下也能做出精确的拟合;其次在本发明技术中,为了克服克里金模型在近似建模中存在部分问题(如相关参数的确定),结合贝叶斯理论,对克里金模型回归系数求解、克里金相关参数求解这两方面进行了研究,进一步改进了克里金模型,增加了克里金模型自适应、自调整能力,降低了拟合结果对数据样本的敏感度,保证了即使有异常点的情况下也能做到拟合结果不受影响。
2)对高程测试点样本进行试验设计优化
在拟合技术中,已有的数据样本的选取对拟合效果的影响非常大。在进行桥面高程拟合时,未见已有技术对高程测试点(测试位置)样本进行优化设计,已有技术一般按分位(如沿跨1/2位置、1/4位置、1/8位置等)布置高程测试点,并不考虑样本分布的均匀性。本发明技术采用先优化设计-再布点的方式。利用数论方法对初始样本点进行试验设计,得到确定性初始样本值,提高了初始试验样本的均匀性和稳健性,在不增加测量工作量的情况下,提高了高程测试点样本质量,对保证后期的模型拟合精度提供了良好基础。
本发明相对于现有技术具有如下的优点:
1、小样本下的高精度拟合。本发明综合利用了克里金模型的技术优势(具有平滑效应以及估计方差最小的统计特征,是一种最优的线性无偏估计),并将其应用到桥面高程拟合中,可实现在较少高程测量数据点的情况下,对全桥面进行全局兼顾局部的高精度拟合。
2、优化了高程测试点样本,进一步提高了拟合精度。本发明不同于传统的按分位布置高程测试点的方式,而是采用基于数论法的样本试验设计对测试点进行了优化,实现了在不增加工作量的情况下提高测试点样本质量,进一步提高模型的预测精度。
3、预测模型具有高自适应性,对数据样本的敏感度低。本发明利用贝叶斯理论对克里金模型回归系数、克里金相关参数进行了求解,改进了克里金模型,建立了贝叶斯-克里金模型,降低了预测模型对模型参数的依赖性,大大提高了模型预测的稳健性。可实现当存在有异常点的情况时,也能保证拟合的精确度。
附图说明
图1为基于贝叶斯理论的克里金模型回归系数求解示意图;
图2为范-克里金模型转换高斯过程示意图;
图3为基于贝叶斯三层先验的克里金回归系数求解示意图;
图4为基于贝叶斯理论的克里金模型相关参数求解示意图;
图5为基于贝叶斯两层先验的克里金相关参数建模示意图;
图6为克里金模型回归系数求解、克里金相关参数求解的综合计算流程图;
图7为GP点集示意图;
图8为本发明的基于贝叶斯-克里金模型的桥面高程拟合方法的流程图;
图9为本发明的基于贝叶斯-克里金模型的桥面高程拟合的效果图;
具体实施方式
下面结合附图和实施例对本发明作进一步说明。
如图1-图9所示的基于贝叶斯-克里金模型的桥面高程拟合方法,包括以下步骤:
S1、建立贝叶斯-克里金拟合模型;
一般来说,克里金模型包含两部分:多项式和随机分布,即:
y(x)=F(β,x)+z(x) (2-1)
式中:
Figure BDA0001829890230000081
其中β为回归系数,f(x)为变量x的多项式函数,p为fi(x)的数目,类似于响应面法中的多项式形式。在设计空间中f(x)提供模拟的全局近似,而z(x)提供模拟的局部近似。z(x)是一个随机过程,服从正态分布N(0,σ2),但是协方差非零,即z(x)不独立,z(x)的协方差矩阵为:
Figure BDA0001829890230000082
式中:
Figure BDA0001829890230000083
为样本点中任何两个样本点w和x的空间相关函数,它对模拟的精确程度起决定性作用,θ为相关函数的参数,式中核函数
Figure BDA0001829890230000084
有多种形式,如高斯函数:
Figure BDA0001829890230000091
克里金模型具有优秀的预测值估计性能及独有的预测均方差(MSE)估计能力,被认为是一种极具潜力的响应面代理模型。但克里金模型在近似建模中仍存在部分问题待完善,如相关参数的确定、优化等。
S11、如图1所示,建立贝叶斯-克里金模型回归系数求解子系统:
惩罚似然是研究拟合近似模型的一种有效方法,例如惩罚最小二乘方法:
Figure BDA0001829890230000092
式中,||y-f(x)Tβ||表示欧几里得范数;γ表示正则化参数或者调谐参数;p(γ,|βi|)是一个事先给定的非负惩罚函数。最小化问题(2-5)中参数β的估计将主要取决于惩罚函数p(γ,|βi|)的选取,不同的惩罚函数得到的解的形式也不同。在下文中,可以看到罚似然方法与贝叶斯理论框架本质是相同的。贝叶斯理论与克里金模型相结合的方法对于求解拟合模型非常有效。在现实中,克里金模型中反应输入之间相关性的相关参数是未知的,并且一般情况下没有途径可以获取它的先验信息。
S111、如图2所示,以范-克里金模型y(x)=f(x)Tβ+z(x)为基础,式中:f(x)=[f1(x)…fp(x)]表示已经选择的基函数,β为回归系数,z(x)是服从均质为0、协方差矩阵为σ2R的高斯过程;得到y*=F*β+z*,此时z*是均值为0、协方差矩阵为σ2I单位矩阵的高斯过程,具体过程如下:
这里提出一种基于多层先验的克里金建模方法,考虑下列模型:
y(x)=f(x)Tβ+z(x) (2-6)
该模型称为范-克里金模型,(当式中f(x)Tβ简化为常数μ时,称为标准-克里金模型)其中,f(x)=[f1(x)…fp(x)]表示已经选择的基函数,βi为回归系数,z(x)是服从均值为零、协方差矩阵为σ2R的高斯分布。由于R是高斯型相关函数,相关矩阵是对称正定的,可利用Cholesky分解将R写为因式分解形式:
R=ΛΛ′ (2-7)
式(2-7)中,Λ是上三角矩阵,称为Cholesky因子。在式(2-6)两边同乘以Λ-1,可得:
Figure BDA0001829890230000101
可知E(Λ-1z)=0,E[(Λ-1z)(Λ-1z)′]=σ2I,定义
Figure BDA0001829890230000102
F*=Λ-1F,z*=Λ-1z,将式(2-8)改写为:
y*=F*β+z* (2-9)
此时z*是均值为0协方差矩阵为σ2I单位矩阵的高斯过程。
由于模型拟合的目的就是寻找一个或几个特定的输入变量对某个感兴趣的输出所发挥的作用以及输入与输出之间具体的关系。因此,变量的选取在模型拟合的分析中是举足轻重的。模型拟合试验的目标是筛选出参数声中的重要变量或者求得模型的稀疏解,由此得到预测可靠且形式简单的拟合模型。在贝叶斯公式里,若对模型中感兴趣的参数赋予一个具有稀疏性的先验,那么这个感兴趣的参数是具有稀疏性的。近年来,Laplace密度函数是一种被经常采用的稀疏性先验,密度函数定义如下:
Figure BDA0001829890230000103
S112、如图3所示,基于贝叶斯三层先验的克里金回归系数建模:
这里将共轭先验稀疏建模方法,用于克里金建模问题。首先,对方差的逆(1/σ2)赋予Gamma先验,即:
p(a|u,v)=Γ(a|u,v) (2-11)
式中:a=1/σ2,u,v是Gamma先验中的超参数。对回归系数β赋予如下变方差的稀疏先验:
Figure BDA0001829890230000104
式中:N(βi|0,τi)是每一个βi对应的一个满足均值为零,方差的为τi的高斯密度函数,而
Figure BDA0001829890230000105
再对参数τi赋予Gamma先验,作为贝叶斯稀疏建模的第二层超先验,即:
p(τi|λ)=Γ(τi|1,λ/2)=0.5λexp(-0.5λτi) (2-13)
因此,在给定参数λ的情况下,可以得到β的先验概率密度:
Figure BDA0001829890230000111
可见,p(β|λ)实际上对应的是具有重尾性的Laplace概率密度函数,能够刻画回归参数β的稀疏性。对参数λ赋予Gamma先验,作为整个稀疏建模过程的第三层超先验,即:
p(λ|η)=Γ(λ|0.5η,0.5η) (2-15)
其中η表示超参数。当给定y*的似然函数p(y*|β,a)、p(β|τ)、p(τi|λ)、p(λ)、p(a)后,通过贝叶斯公式不难得到β的后验分布:
p(β|y*,τ,a,λ)∝p(y*|β,a)p(β|τ)p(τ|λ)p(λ)p(a)∝p(y*|β,a)p(β|λ)p(λ)。
可以得到β的最大后验概率估计为:
Figure BDA0001829890230000112
经过计算,可以发现p(β|y*,τ,a,λ)对应的是一个多维正态分布,即:
p(β|y*,τ,a,λ)~N(β|[aF*′F*+Ω(τ)]-1aF*′y*,[aF*′F*+Ω(τ)]-1) (2-17)
根据最大后验概率准则,参数β的估计值为:
Figure BDA0001829890230000113
基于贝叶斯理论的克里金模型回归系数求解过程如图1所示。
S12、建立贝叶斯-克里金模型相关参数求解子系统:
S121、如图4所示,以范-克里金模型y(x)=f(x)Tβ+z(x)为基础,输出y的似然函数为:
Figure BDA0001829890230000114
S122、如图5所示,基于贝叶斯两层先验的克里金相关参数建模:
本节贝叶斯建模的重点在于参数θ。这里提出了关于参数θ的一种基于两层先验的克里金建模方法,即:第一层为相关参数θ赋予变方差的高斯分布,第二层为变方差赋予Jeffreys非信息超先验。首先对回归系数β和方差逆a赋予如下常数先验:
Figure BDA0001829890230000121
为了得到自适应解,对相关参数θ提出一个多层先验:
Figure BDA0001829890230000122
式中,N(θi|0,ωi)是一个变方差的高斯分布,ωi表示方差,
Figure BDA0001829890230000123
然后对每一个ωi取非信息的Jeffreys超先验:p(ωi)=1/ωi。非信息Jeffreys超先验是不含有参数的先验,并且是非正常先验。因此,由贝叶斯多层先验方法得到的后验分布是:
Figure BDA0001829890230000124
对(2-22)式取自然对数,可以得到:
Figure BDA0001829890230000125
由于输出y的罚似然函数是:
Figure BDA0001829890230000126
由式(2-23)可知,θ,a的最大后验概率估计可最终通过最大化下式得到:
Figure BDA0001829890230000127
由于后验分布不是标准形式,因此给实际计算带来困难,虽然目标是求出θ,a的最大后验估计,但是由于变方差ω=(ω1,…ωL)也是未知的,因此在θ,a的同时,还需要估计ω=(ω1,…ωL),这里将ω=(ω1,…ωL)看作是缺失数据,利用EM算法进行计算,限于篇幅,计算过程不再赘述,其中EM算法为本领域常见的算法之一。
综合克里金模型回归系数求解、克里金相关参数求解的整个综合计算流程如图6所示。
S2、高程测试点样本试验优化;
初始训练样本点(高程测试点样本)的选取对代理模型的预估精度影响较大,为了提高模拟精度,这里摒弃以前传统的按分位(如沿跨1/2位置、1/4位置、1/8位置等)布置高程测试点位置,对高程测试点采集位置进行试验设计优化。常用的试验设计方法主要有正交实验设计,均匀设计,中心复合设计,拉丁超立方抽样等,本发明采用数论方法对高程测试点进行试验设计优化。数论方法的主要特点就体现在样本点的产生上,它提出了一系列产生样本点的确定性方法,特点是,舍弃样本的随机性,注重样本的均匀性。它的主要优点有:样本的均匀性好,计算效率高且精度有保证;由于样本的产生是确定性方法,计算结果也是确定值,因此具有较好的稳健性。
数论方法采用偏差为测度衡量一组样本点的均匀性,令s为一正整数及Cs表示s维空间的单位立方体,0≤x1≤1,...,0≤xs≤1,令n为一正整数及
Figure BDA0001829890230000131
表示Cs中的一个点集,其中
Figure BDA0001829890230000132
对于任意δ=(δ1,...,δs)∈Cs,令Nn(δ,Jn)代表满足
Figure BDA0001829890230000133
条件的点的数量,则:
Figure BDA0001829890230000134
则称点集Jn有偏差D(n,Jn)。如果:D(n,Jn)=O(n-1/2) (2-27)
成立,则认为Jn是均匀分布于Cs的数论点集。有许多种方法可以得到数论点集,例如GLP点集、Halton点集、Halton Leaped点集、Hamersley点集、GP点集。本文采用GP点集,其最大优点是对高维问题具有很好的适用性。由一个所谓的好点得到的集合称为GP(GoodPoint)点集,在使用时一般采用平方根序列、分圆域方法等产生好点的方法。
令δ={δ1,...,δs}∈Cs,如果点集:Jn(k)={[δ1k],...,[δsk]},1≤k≤n (2-28)
具有偏差:D(n,Jn)≤O(δ,ε)n-1+ε (2-29)
则称Jn点集为一个GP点集,而δ为一个好点,其中[.]表示取小数运算符。在实际使用中,常用的是分圆域方法来得到δ:
Figure BDA0001829890230000135
其中h为素数且h≥2s+3,S为维数,该集合具有的偏差为:
Figure BDA0001829890230000136
图7为二维变量设计样本为20组的GP点集示意图,从图中可以看到数论点集具有较好的均匀性。GP点集的生成机制是确定的,具有稳健性,操作也是十分方便的,点集的维数和容量都可以方便的扩张。
在高程测试点样本采集阶段,对测试点样本进行优化设计,采用确定性的数论法进行初始训练样本试验设计,保证了计算方法的稳健性,而且数论抽样方法产生的样本具有优秀的均匀性。利用样本的空间均匀性特点,可即有效减少测试工作量,又能提高样本质量。按照数论法优化设计得到的高程测试点位置,在桥梁现场进行高程测试。
S3、建立贝叶斯-克里金预测模型;
利用采集得到的高程测试点样本对贝叶斯-克里金拟合模型进行训练,得到最终的贝叶斯-克里金预测模型。
S4、通过S3中得到的贝叶斯-克里金预测模型对全桥面进行高程拟合评估。
利用建立好的贝叶斯-克里金预测模型,输入其他桥面点平面位置(xy),拟合估算出对应的高程,并利用克里金预测模型预测出均方差对拟合效果进行评估。本发明的总流程图参见8,经本发明的桥面高程拟合效果图参见图9。
本发明将贝叶斯和克里金结合:基于非信息超先验,对克里金模型的基函数系数、相关参数赋予了多层先验约束,利用EM算法,求解基函数系数、相关参数的最大后验估计,对克里金模型进行改进,建立贝叶斯-克里金模型,增强了模型的自适应性和稳健性。同时进行基于数论法的高程测量样本的试验优化设计。最后将其成功的应用到了桥面高程拟合中,实现了可依靠较少数量的已有高程测试点数据,以较小的计算量,对全桥面实现高精度拟合。
上述具体实施方式为本发明的优选实施例,并不能对本发明进行限定,其他的任何未背离本发明的技术方案而所做的改变或其它等效的置换方式,都包含在本发明的保护范围之内。

Claims (3)

1.一种基于贝叶斯-克里金模型的桥面高程拟合方法,其特征在于,包括以下步骤:
S1、建立贝叶斯-克里金拟合模型;
S11、建立贝叶斯-克里金模型回归系数求解子系统:
S111、以范-克里金模型y(x)=f(x)Tβ+z(x)为基础,式中:f(x)=[f1(x)L fp(x)]表示己经选择的基函数,β为回归系数,z(x)是服从均质为0、协方差矩阵为σ2R的高斯过程;得到y*=F*β+z*,此时z*是均值为0、协方差矩阵为σ2I单位矩阵的高斯过程;
S112、基于贝叶斯三层先验的克里金回归系数建模:
第一层为对随机过程的方差逆(1/σ2)赋予Gamma先验:p(a|u,v)=Γ(a|u,v),式中:a=1/σ2,u,v是Gamma先验中的超参数;
对回归系数β赋予变方差的高斯先验:
Figure FDA0002794528780000011
式中:N(βi|0,τi)是每一个βi对应的一个满足均值为零,参数τi为高斯密度函数,
Figure FDA0002794528780000012
第二层为对参数τi赋予Gamma先验:
p(τi|λ)=Γ(τi|1,λ/2)=0.5λexp(-0.5λτi),
第三层为对参数λ赋予Gamma先验:
p(λ|η)=Γ(λ|0.5η,0.5η),式中:η表示超参数;
得出β的估值为:
Figure FDA0002794528780000021
S12、建立贝叶斯-克里金模型相关参数求解子系统:
S121、以范-克里金模型y(x)=f(x)Tβ+z(x)为基础,输出y的似然函数为:
Figure FDA0002794528780000022
S122、基于贝叶斯两层先验的克里金相关参数建模:
第一层为相关参数θ赋予变方差的高斯分布;
首先对回归系数β和方差逆a赋予如下常数先验:
Figure FDA0002794528780000023
为了得到自适应解,对相关参数θ提出一个多层先验:
Figure FDA0002794528780000024
式中:N(θi|0,ωi)是一个变方差的高斯分布,ωi表示方差,
Figure FDA0002794528780000025
第二层为对每一个ωi取Jeffreys非信息的超先验:
p(ωi)=1/ωi
S2、高程测试点样本试验优化:采用数论方法对高程测试点样本试验进行优化;
S3、建立贝叶斯-克里金预测模型;
S4、对全桥面进行高程拟合评估。
2.根据权利要求1中所述的基于贝叶斯-克里金模型的桥面高程拟合方法,其特征在于:所述S4包括以下内容:
向S3的贝叶斯-克里金预测模型中,输入其他桥面点的平面位置,拟合估算出对应的高程。
3.根据权利要求2中所述的基于贝叶斯-克里金模型的桥面高程拟合方法,其特征在于:所述S4还包括利用S3的贝叶斯-克里金预测模型预测均方差对拟合效果进行评估。
CN201811200557.4A 2018-10-16 2018-10-16 一种基于贝叶斯-克里金模型的桥面高程拟合方法 Active CN109540089B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811200557.4A CN109540089B (zh) 2018-10-16 2018-10-16 一种基于贝叶斯-克里金模型的桥面高程拟合方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811200557.4A CN109540089B (zh) 2018-10-16 2018-10-16 一种基于贝叶斯-克里金模型的桥面高程拟合方法

Publications (2)

Publication Number Publication Date
CN109540089A CN109540089A (zh) 2019-03-29
CN109540089B true CN109540089B (zh) 2021-05-14

Family

ID=65843751

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811200557.4A Active CN109540089B (zh) 2018-10-16 2018-10-16 一种基于贝叶斯-克里金模型的桥面高程拟合方法

Country Status (1)

Country Link
CN (1) CN109540089B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112598178B (zh) * 2020-12-23 2022-03-29 北方民族大学 一种基于贝叶斯框架下的大气污染物浓度预测方法
CN113094843B (zh) * 2021-04-30 2022-08-09 哈尔滨工业大学 一种基于贝叶斯网络的梁式桥评估的条件概率的求解方法
CN115017767A (zh) * 2022-06-02 2022-09-06 厦门大学 基于贝叶斯正则化的桥梁影响线识别与不确定性量化方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7961922B2 (en) * 2007-05-31 2011-06-14 The Board Of Regents Of The University Of Texas System Systems and methods for processing medical image data to facilitate comparisons among groups of subjects
CN103530904A (zh) * 2013-11-04 2014-01-22 东南大学 一种基于克里金方法的水下地形数字高程建立方法
CN103617462A (zh) * 2013-12-10 2014-03-05 武汉大学 一种基于地理统计学的风电场风速时空数据建模方法
CN104316978A (zh) * 2014-10-29 2015-01-28 中国石油天然气股份有限公司 地球物理的近地表三维速度场研究方法和装置
CN107423561A (zh) * 2017-07-11 2017-12-01 云南瀚哲科技有限公司 一种土壤属性插值的估算方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7961922B2 (en) * 2007-05-31 2011-06-14 The Board Of Regents Of The University Of Texas System Systems and methods for processing medical image data to facilitate comparisons among groups of subjects
CN103530904A (zh) * 2013-11-04 2014-01-22 东南大学 一种基于克里金方法的水下地形数字高程建立方法
CN103617462A (zh) * 2013-12-10 2014-03-05 武汉大学 一种基于地理统计学的风电场风速时空数据建模方法
CN104316978A (zh) * 2014-10-29 2015-01-28 中国石油天然气股份有限公司 地球物理的近地表三维速度场研究方法和装置
CN107423561A (zh) * 2017-07-11 2017-12-01 云南瀚哲科技有限公司 一种土壤属性插值的估算方法

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
EFFICIENT OPTIMIZATION DESIGN METHOD USING KRIGING MODEL;Shinkyu Jeong;《Journal of Aircraft》;20050930;1-11 *
Why do we need and how should we implement Bayesian kriging methods;Jurgen Pilz;《Stoch Environ Res Risk Assess》;20070620;621-632 *
克里金模型及其在全局优化设计中的应用;于向军;《中国工程机械学报》;20060731;第4卷(第3期);259-261 *
基于动态贝叶斯模型的大悬臂展翅宽箱梁桥性能预测评估;贾布裕;《铁道科学与工程学报》;20160731;第13卷(第7期);1303-1308 *
基于贝叶斯克里金的地下空间多源数据建模;李晓军;《同济大学学报(自然科学版)》;20140331;第42卷(第3期);406-412 *

Also Published As

Publication number Publication date
CN109540089A (zh) 2019-03-29

Similar Documents

Publication Publication Date Title
Sirignano et al. DGM: A deep learning algorithm for solving partial differential equations
Durrant winGamma: A non-linear data analysis and modelling tool with applications to flood prediction
CN110909926A (zh) 基于tcn-lstm的太阳能光伏发电预测方法
CN111860982A (zh) 一种基于vmd-fcm-gru的风电场短期风电功率预测方法
CN109540089B (zh) 一种基于贝叶斯-克里金模型的桥面高程拟合方法
CN108304679A (zh) 一种自适应可靠性分析方法
Jayaram et al. A Neural Network Approach To Tolerance Synthesis and Cost Optimisation in Assembly
CN112418482A (zh) 一种基于时间序列聚类的云计算能耗预测方法
CN115203865B (zh) 一种基于数字孪生的产品装配过程机械性能在线预测方法
CN112818595A (zh) 一种火电厂蒸发区的数字孪生模型数据的修正方法及系统
CN114707712A (zh) 一种发电机组备件需求的预测方法
US20230062600A1 (en) Adaptive design and optimization using physics-informed neural networks
CN112257021A (zh) 一种可选的克里金空间插值降雨量估算方法
US8813009B1 (en) Computing device mismatch variation contributions
Mattey et al. A physics informed neural network for time-dependent nonlinear and higher order partial differential equations
CN115983029A (zh) 一种航空装备可靠性仿真数字孪生模型构建方法、设备、介质
Regazzoni et al. A physics-informed multi-fidelity approach for the estimation of differential equations parameters in low-data or large-noise regimes
CN116432543A (zh) 功率半导体模块剩余寿命预测方法、终端设备及存储介质
Drakoulas et al. FastSVD-ML–ROM: A reduced-order modeling framework based on machine learning for real-time applications
CN115421216A (zh) 一种基于stl-arima-nar混合模型的中长期月降雨预报方法
CN108228978A (zh) 结合互补集合经验模态分解的Xgboost时间序列预测方法
Wang et al. Adaptive echo state network with a recursive inverse-free weight update algorithm
CN112381279B (zh) 一种基于vmd和bls组合模型的风电功率预测方法
CN112232565A (zh) 基于两阶段的时间序列预测方法、预测系统、终端及介质
CN116776695A (zh) 基于超高阶有限元技术的一维电磁计算方法、系统及设备

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant