对大样本空间复杂概率分布的无自相关采样方法技术

本发明专利技术公布了一种对大样本空间复杂概率分布的无自相关采样方法，目的是解决MCMC方法的自相关问题以及跳步采样的效率损失问题。技术方案是：在产生样本时，通过设置样本缓冲区，消除样本之间的关联；或者作为后处理方法，按顺序读取已产生的样本，通过样本缓冲区更新样本序列，以消除样本序列自相关。每次产生一个样本，或者从文件中读取到一个样本时，将得到的样本存储到样本缓冲区中，并在缓冲区充满时从缓冲区中随机选择样本进行输出。通过这种随机输出的方式在样本之间增加了一层额外的随机性，最终达到消除样本序列自相关的目的。采用本发明专利技术计算速度快，维持了MCMC效率，额外增加的缓冲区操作的开销相对可忽略，同时消除了样本自相关。

A Non-autocorrelation Sampling Method for Complex Probability Distribution in Large Sample Space

【技术实现步骤摘要】
对大样本空间复杂概率分布的无自相关采样方法
本专利技术涉及对大样本空间复杂概率分布样本的采样方法，尤指基于马尔科夫链的无自相关采样方法。
技术介绍
计算机对随机事件的处理是工程应用分析的必要手段。近几十年来，随着计算机性能的飞速发展，计算机已广泛应用于航空航天、汽车船舶设计制造、桥梁建筑设计制造、环境工程、气象预报、高分子材料等方面，而在这些工程应用中，均需要通过计算机对随机事件进行计算分析，以达到如桥梁可靠性检测、气象预报等目的。在随机事件分析中，一个必要步骤是对复杂概率分布进行采样，以获取独立、随机的样本。高效的采样算法是提高采样效率的必要手段。一般而言，工程应用中涉及到的随机事件来源复杂，数据量大，样本概率计算难度大，样本种类数可能多达1020量级以上。此类随机事件的分析需要大量的计算资源。除了尽可能提高硬件性能以满足计算资源的需求外，对采样方法的优化必不可少。大样本空间复杂概率分布采样问题可定义如下：·样本空间X＝{x1,x2,…,xi,…,xN}表示可能存在的N个样本(1≤i≤N)。N的值十分巨大，可能高达1020量级。·复杂概率分布由与领域相关的复杂概率分布函数f(x)描述，即采样得到样本xi的概率为f(xi)。概率值f(xi)的计算是困难的(一般为指数时间复杂度)，需要大量的计算资源。·构造样本序列t1,t2,…,tj,….,tM，其中对于任意的j(1≤j≤M，M的值随应用的不同而不同，1≤M≤N)均有样本tj∈X，并使得序列中的样本服从概率分布函数f描述的概率分布。概率分布函数f的计算一般需要大量的计算资源，是采样过程的关键点和难点。目前常用的采样方法包括暴力采样方法、拒绝采样方法(A.Neville等，2017年10月，《自然-物理》,具有超越近期物理实验性能的经典玻色采样算法)及马尔科夫链蒙特卡洛方法(J.S.Liu，1996年6月，《统计与计算》,Metropolized独立采样方法以及与拒绝采样方法和重要性采样方法的对比)。1.暴力采样方法暴力采样方法通过先完整地计算整个待采样的概率分布，而后完成采样。而计算整个待采样概率分布的计算资源需求十分巨大，不适用于工程应用中。2.拒绝采样方法拒绝采样方法首先基于一个假定的简单的辅助概率g分布来生成待定样本(记为x)，再计算待定样本x在实际的待采样概率分布f中的概率值f(x)。辅助概率g需要乘以一个因子λ使得对于任意样本均有λ·g(x)>f(x)。此后，待定样本的实际概率与辅助概率的比值f(x)/λ·g(x)将作为将待定样本进行输出的依据，即以f(x)/λ·g(x)为概率进行输出，或以1–f(x)/λ·g(x)为概率将待定样本舍弃(即拒绝)。虽然拒绝采样方法在产生样本时不需要计算整个概率分布，但仍然有以下缺点：(1)由于实际概率分布的概率值计算困难，其最大概率值获取困难，因此因子λ的值难以确定；(2)当辅助概率分布定位过小时，会对实际概率分布有一定的截断效果，使得采样结果失真；(3)当辅助概率分布定位过大时，接受概率很低，样本被拒绝的概率大大提升，使得采样效率大幅降低；(4)即便较好地确定了因子λ的值，平均每接受一个样本需要计算的概率数随待采样概率分布的复杂度提升而提高，对计算资源的需求依然较大，效率不高。3.马尔科夫链蒙特卡洛(MCMC)MCMC方法通过构造一条以待采样概率分布为收敛概率分布的马尔科夫链来进行采样。每计算一个概率值时，马尔科夫链便能扩展一个节点，同时可进行一个样本输出，是目前最高效的采样方法。马尔科夫链蒙特卡洛方法的总体流程如图1所示，包括。第一步：选取初态。1.1定义马尔科夫链的状态空间为S＝{s1,s2,…,si,…,sN}，其中对于正整数i(1≤i≤N)，有si与样本空间中X＝{x1,x2,…,xi,…,xN}的元素xi对应；1.2构造基于马尔科夫链状态空间S上的样本概率函数fs(s)，使得对于整数i，均有fs(si)＝f(xi)，其中f(xi)表示采样得到样本xi的概率。1.3任意选取状态空间S上的易采样随机分布g(s)作为辅助概率分布(例如均匀分布等)，其中g(si)表示采样得到状态si作为样本的概率。从g(s)的采样是简单的，一般选取独立采样方法，即从g(s)中的采样过程与当前状态无关；1.4按均匀随机分布产生正整数o∈[1,N]，并以so作为马尔科夫链的初态，计算状态so对应的样本概率po＝fs(so)，并输出状态so的对应样本xo；1.5设置马尔科夫链当前状态sc＝so；1.6根据用户需求与应用特点，设置采样样本数量目标值Ns。第二步：产生待选状态2.1从条件概率分布g(s)中随机选取状态sn(1≤n≤N)作为待选状态；2.2计算待选状态对应的概率值fs(sn)。第三步：决策确定后续状态，输出样本3.1计算接受概率paccept＝min(1,(fs(sn)g(sc))/(fs(sc)g(sn)))；3.2产生[0,1]区间内的随机数u；3.3若u<paccept，令sc＝sn，即接受待选状态sn为马尔科夫链的新当前状态；否则sc保持不变，即拒绝待选状态；3.4生成并输出状态sc对应的样本。第四步：判断采样是否完成若采样得到的样本数已经达到既定目标值，跳转第五步；若尚未达到既定目标值，跳转到到第二步。第五步：结束。MCMC方法能够在每计算一个样本值时即能产生一个样本，然而产生的样本序列存在一个很大的缺点，即由于第二步的接受或拒绝过程的事件概率严重依赖于当前样本，使得产生的序列中可能会出现大量的连续重复样本，样本序列中存在严重的自相关。为了定量表征样本序列的自相关程度，一般通过各阶自相关系数，尤其是一阶自相关系数r1来评估样本序列的自相关程度，如式(1)所示。其中k表示评估的自相关阶数，xi为第i个样本值，为样本均值。设置k＝1则对一阶自相关进行评估。r1的值介于-1与1之间，其绝对值越接近1，意味着样本序列中的自相关程度越大，r1值为0则意味着样本之间相互独立，即无自相关。而通过马尔科夫链产生的样本序列其一阶自相关系数r1可能会高达0.8以上。目前常用的用于消除自相关的方法称为“跳步采样”，是基本MCMC方法的一个改进，其总体流程如图2所示，包括：第一步：选取初态1.1定义马尔科夫链的状态空间为S＝{s1,s2,…,si,…,sN}，其中对于正整数i(1≤i≤N)，有si与样本空间中X＝{x1,x2,…,xi,…,xN}的元素xi对应；1.2构造基于马尔科夫链状态空间S上的样本概率函数fs(s)，使得对于整数i，均有fs(si)＝f(xi)，其中f(xi)表示采样得到样本xi的概率。1.3任意选取状态空间S上的易采样随机分布g(s)作为辅助概率分布(例如均匀分布等)，其中g(si)表示采样得到状态si作为样本的概率。从g(s)的采样是简单的，一般选取独立采样方法，即从g(s)中的采样过程与当前状态无关；1.4按均匀随机分布产生正整数o∈[1,N]，并以so作为马尔科夫链的初态，计算状态so对应的样本概率po＝fs(so)，并输出状态so的对应样本xo；1.5设置马尔科夫链当前状态sc＝so；1.6确定跳步采样的步长k。一般k取值为100；1.7根据用户需求与应用特点，设置采样样本本文档来自技高网...

【技术保护点】
1.一种对大样本空间复杂概率分布的无自相关采样方法，其特征在于包括以下步骤：第一步：选取初态，初始化缓冲区，方法是：1.1定义马尔科夫链的状态空间为S＝{s1,s2,…,si,…,sN}，其中对于正整数i，1≤i≤N，有si与样本空间中X＝{x1,x2,…,xi,…,xN}的元素xi对应；1.2构造基于马尔科夫链状态空间S上的样本概率函数fs(s)，使得对于整数i，均有fs(si)＝f(xi)，其中f(xi)是复杂概率分布函数，表示采样得到样本xi的概率；1.3任意选取状态空间S上的易采样随机分布g(s)作为辅助概率分布，其中g(si)表示采样得到状态si作为样本的概率；1.4设置样本缓冲区容量L初始值，L为正整数，初始化样本缓冲区为空，容量为可容纳L个样本，其中每个可存储样本的位置称为一个缓冲槽，初始化当前存储位置变量p＝1；1.5按均匀随机分布产生正整数o∈[1,N]，并以so作为马尔科夫链的初态，计算状态so对应的样本概率po＝fs(so)，并输出状态so的对应样本xo；1.6设置马尔科夫链当前状态sc＝so；1.7根据用户需求与应用特点，设置采样样本数量目标值Ns；第二步：产生待选状态，方法是：2.1从g(s)中随机采样得到状态sn作为待选状态，1≤n≤N；2.2计算待选状态对应的概率值fs(sn)；第三步：决策确定后续状态，生成样本，方法是：3.1计算接受概率paccept＝min(1,(fs(sn)g(sc))/(fs(sc)g(sn)))；3.2产生[0,1]区间内的随机数u；3.3若u<paccept，则令sc＝sn；否则sc保持不变；3.4生成状态sc对应的样本，并记录该样本产生的序号ngenerate；第四步：判定使用方式，方法是：若对样本实时性要求高，且无充分时间进行调试，设置模式变量mode＝1，并转第八步；若所进行的是调试过程，即测量所需缓冲区容量大小，则设置模式变量mode＝2，转第五步；第五步：存储样本：将第3.4步产生的样本存储到样本文件中；第六步：判断采样是否完成：若第3.4步得到的样本产生的序号ngenerate＝Ns，转第七步；若第3.4步得到的样本产生序号ngenerate<Ns，则转第二步；第七步：读取样本：从样本文件中按顺序读取一个样本；第八步：判断样本缓冲区是否充满：检查样本缓冲区最后一个缓冲槽是否有样本，若有，转第九步；若没有，转第十步；第九步：随机输出样本，方法是：9.1产生[1,L]之间的随机正整数rb；9.2将样本缓冲区中的第rb个样本缓冲槽中存储的样本输出，并清空第rb个缓冲槽，记录样本的输出序号noutput；9.3若mode＝1，将第3.4步中生成的样本存储到样本缓冲区的第rb个缓冲槽中；若mode＝2，则将第七步中读取的样本存储到样本缓冲区的第rb个缓冲槽中；9.4转第十一步；第十步：填充样本缓冲区10.1若mode＝1，则将第3.4步中产生的样本存储到第p个样本缓冲槽中；10.2若mode＝2，则将第七步中读取的样本存储到第p个样本缓冲槽中；10.3p＝p+1；10.4转第十一步；第十一步：判定使用方式：若mode＝1，转第十二步；若mode＝2，则转第十三步；第十二步：判断采样是否完成：若第9.2步记录的输出样本序号noutput＝Ns–L，转第十四步；若第9.2步记录的输出样本序号noutput<Ns–L，则转第二步；第十三步：判断样本是否读取完成：若样本文件中的样本已完全被读出，转第十四步；若文件中仍然有未读出的样本，则转第七步；第十四步：清空样本缓冲区，方法是：14.1初始化剩余样本数变量Nleft＝L；14.2产生[1,Nleft]之间的随机正整数rc，并将缓冲区中第rc个缓冲槽中样本输出；14.3将Nleft个缓冲槽中样本移动到第rc个缓冲槽；14.4Nleft＝Nleft–1；14.5若Nleft值为0，转第十五步，否则转14.2；第十五步：结束。...

【技术特征摘要】
1.一种对大样本空间复杂概率分布的无自相关采样方法，其特征在于包括以下步骤：第一步：选取初态，初始化缓冲区，方法是：1.1定义马尔科夫链的状态空间为S＝{s1,s2,…,si,…,sN}，其中对于正整数i，1≤i≤N，有si与样本空间中X＝{x1,x2,…,xi,…,xN}的元素xi对应；1.2构造基于马尔科夫链状态空间S上的样本概率函数fs(s)，使得对于整数i，均有fs(si)＝f(xi)，其中f(xi)是复杂概率分布函数，表示采样得到样本xi的概率；1.3任意选取状态空间S上的易采样随机分布g(s)作为辅助概率分布，其中g(si)表示采样得到状态si作为样本的概率；1.4设置样本缓冲区容量L初始值，L为正整数，初始化样本缓冲区为空，容量为可容纳L个样本，其中每个可存储样本的位置称为一个缓冲槽，初始化当前存储位置变量p＝1；1.5按均匀随机分布产生正整数o∈[1,N]，并以so作为马尔科夫链的初态，计算状态so对应的样本概率po＝fs(so)，并输出状态so的对应样本xo；1.6设置马尔科夫链当前状态sc＝so；1.7根据用户需求与应用特点，设置采样样本数量目标值Ns；第二步：产生待选状态，方法是：2.1从g(s)中随机采样得到状态sn作为待选状态，1≤n≤N；2.2计算待选状态对应的概率值fs(sn)；第三步：决策确定后续状态，生成样本，方法是：3.1计算接受概率paccept＝min(1,(fs(sn)g(sc))/(fs(sc)g(sn)))；3.2产生[0,1]区间内的随机数u；3.3若u<paccept，则令sc＝sn；否则sc保持不变；3.4生成状态sc对应的样本，并记录该样本产生的序号ngenerate；第四步：判定使用方式，方法是：若对样本实时性要求高，且无充分时间进行调试，设置模式变量mode＝1，并转第八步；若所进行的是调试过程，即测量所需缓冲区容量大小，则设置模式变量mode＝2，转第五步；第五步：存储样本：将第3.4步产生的样本存储到样本文件中；第六步：判断采样是否完成：若第3.4步得到的样本产生的序号ngenerate＝Ns，转第七步；若第3.4步得到的样本产生序号ngenerate<Ns，则转第二步；第七步：读取样...

【专利技术属性】
技术研发人员：吴俊杰，刘雍，熊敏，徐平，强晓刚，黄安琪，付祥，邓明堂，
申请(专利权)人：中国人民解放军国防科技大学，
类型：发明
国别省市：湖南,43