[go: up one dir, main page]

CN117574790B - Design method of cross-dimension Bayesian sampler based on physical space tree structure - Google Patents

Design method of cross-dimension Bayesian sampler based on physical space tree structure Download PDF

Info

Publication number
CN117574790B
CN117574790B CN202410078165.4A CN202410078165A CN117574790B CN 117574790 B CN117574790 B CN 117574790B CN 202410078165 A CN202410078165 A CN 202410078165A CN 117574790 B CN117574790 B CN 117574790B
Authority
CN
China
Prior art keywords
tree structure
model
node
nodes
distribution
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
CN202410078165.4A
Other languages
Chinese (zh)
Other versions
CN117574790A (en
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.)
Central South University
Original Assignee
Central South University
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 Central South University filed Critical Central South University
Priority to CN202410078165.4A priority Critical patent/CN117574790B/en
Publication of CN117574790A publication Critical patent/CN117574790A/en
Application granted granted Critical
Publication of CN117574790B publication Critical patent/CN117574790B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/38—Processing data, e.g. for analysis, for interpretation, for correction
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N5/00—Computing arrangements using knowledge-based models
    • G06N5/01—Dynamic search techniques; Heuristics; Dynamic trees; Branch-and-bound
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N7/00—Computing arrangements based on specific mathematical models
    • G06N7/01—Probabilistic graphical models, e.g. probabilistic networks
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/08—Probabilistic or stochastic CAD
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30—Assessment of water resources

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Evolutionary Computation (AREA)
  • Artificial Intelligence (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Mathematical Physics (AREA)
  • Computing Systems (AREA)
  • Computer Hardware Design (AREA)
  • Pure & Applied Mathematics (AREA)
  • Geometry (AREA)
  • Medical Informatics (AREA)
  • Probability & Statistics with Applications (AREA)
  • Algebra (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Computational Linguistics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • Remote Sensing (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • Geophysics (AREA)
  • Complex Calculations (AREA)

Abstract

The invention discloses a design method of a cross-dimensional Bayesian sampler based on a physical space tree structure, which relates to the technical field of geophysical inversion and comprises the following steps: measuring the response of the underground electromagnetic field of a certain area; constructing a tree structure initial model; optional tree structurem t Tree structure change is carried out on the class of nodes of the model to form a new candidate modelm'The method comprises the steps of carrying out a first treatment on the surface of the New candidate modelm'Inversely transforming to obtain a proposal model, and calculating forward modeling response; calculating a likelihood function; calculating acceptance rate of suggestion modelαJudgingαWhether or not it is smaller than a random number, if so, then there ism t+1 = m'If otherwise there ism t+1 =m t The method comprises the steps of carrying out a first treatment on the surface of the Taking outt=t+1Repeated sampling, per runZSubsampling holds the one-time sampling model untilt=QEnding the time; each will beZAnd carrying out statistical calculation on the sampling model stored by sub-sampling to obtain one-dimensional edge probability density distribution. The invention solves the problem that the prior distribution of the nodes of the wavelet domain tree structure is difficult to select, and improves the sampling efficiency.

Description

Design method of cross-dimension Bayesian sampler based on physical space tree structure
Technical Field
The invention relates to the technical field of geophysical inversion, in particular to a design method of a cross-dimensional Bayesian sampler based on a physical space tree structure.
Background
Geophysical magnetotelluric exploration is a geophysical method for acquiring underground physical response by using a natural electromagnetic field source, and is widely applied to research in the fields of oil gas, mineral products, underground deep-layer ground electrical structures and the like. The inversion results also show non-uniqueness due to interference in the actual observations. Therefore, along with the deep exploration and the rapid improvement of the development level of computer software and hardware, the reliability evaluation of the magnetotelluric inversion result gradually becomes the focus and hot spot of the study of domestic and foreign students. Uncertainty analysis is carried out on the inversion result, the distribution situation of possible underground electrical property is given, quantitative evaluation of the inversion result by people is greatly improved, and guidance is given in actual engineering.
At present, probability inversion based on a Bayesian framework is performed in order to embody uncertainty of an inversion result and quantitatively evaluate the inversion result. Model extraction is carried out by using a Markov Monte Carlo (Markov chain MonteCarlo, mcMC) sampling method, and statistical inference of parameters is carried out.
In bayesian inversion, the dimensions of model parameters cannot be determined in advance because the knowledge of the subsurface structure for actual exploration is limited. Literature (Green, P.J.,1995.Reversible jump Markov chain Monte Carlo computation and Bayesian model determination,Biometrika,82 (4), 711-732) considers model parameters as unknowns and incorporates them into a bayesian framework, a method called trans-dimensional bayesian inversion. The complexity of the model is driven by data through the cross-dimensional Bayesian inversion, and the optimization parameters do not need to be selected manually, so that the influence of human factors on the inversion result is avoided, and the method becomes a mainstream method. In the literature (Mallingerno, A.,2002.Parsimonious Bayesian Markov chain Monte Carlo inversion in a nonlinear geophysical problem,Geophys.J.Int, 151 (3), 675-688) which is a direct current resistance inversion for the first time, the method gives a natural extension of Bayesian parameter estimation, i.e. parameterizing the earth model into a layered medium, inverting the number of layers as parameters, and obtaining samples of the layered medium from posterior distribution. By means of this automatic parameterization, it is possible to provide a higher model resolution in areas with better data coverage and low noise. However, this parameterization method cannot address the high-dimensional magnetotelluric inversion problem.
Currently, the Voronoi cell method is a common method among methods for performing parameter cross-dimensional sampling in a high-dimensional parameter space. Literature (Bodin, t., sambridge, m., rawlinson, N. & aroucau, p.,2012.Transdimensional tomography with unknown data noise,Geophys.J.Int, 189, 1536-1556.) discloses the discretization of subsurface conductivity profiles into a series of Voronoi polygons by forming a set of continuous polygons from perpendicular bisectors connecting adjacent nuclei. The cross-dimensional sampling of the number of model parameters is achieved by increasing or decreasing the number of nuclei. The method has the advantage that the parameter position and physical parameters of the geophysical model space can be represented by specifying the position of atomic nuclei and the attribute value inside each cell, so that the model is flexibly divided. However, the Voronoi cell method requires three unknowns to be contained in one cell. As the data volume of the geophysical model increases, the complexity and time cost of the geophysical model increases. Furthermore, the mesh subdivision of the Voronoi cell subdivision method is a non-uniform arbitrary subdivision, which presents challenges for the finite difference-based geoelectromagnetic forward computing method. Literature (Hawkins, R. & Sambridge, m.,2015.Geophysical imaging using transdimensional trees,Geophys.J.Int,203,972-1000.) proposes a method for sampling a wavelet tree structure. The method converts inversion of the geophysical model in physical space into research and operation in a tree structure based on wavelet basis functions, and can represent a complex geophysical model with a small number of parameters. However, since the wavelet coefficients are stored in the tree-like structure nodes and do not represent any physical significance, the method is difficult to choose prior distribution and cannot correspond to the conductivity range in the physical space. Meanwhile, the selection of different wavelet basis functions can directly influence the final inversion result, so that the application condition of the method is limited to a certain extent. Therefore, it is necessary to provide a new tree structure to solve the above technical problems.
Disclosure of Invention
The invention mainly aims to provide a design method of a cross-dimensional Bayesian sampler based on a physical space tree structure, which aims to solve the problems that prior distribution of nodes of a wavelet domain tree structure is difficult to select and sampling efficiency is low.
In order to achieve the above object, the present invention provides a design method of a cross-dimensional bayesian sampler based on a physical space tree structure, which includes:
s1: measuring by using an electromagnetic measuring instrument to obtain the response of an underground electromagnetic field of a certain area, and determining the measured underground depth;
s2: constructing a tree structure initial model, and constructing a Markov chain, wherein the method specifically comprises the following steps of:
calculating the maximum depth of the tree structure, wherein the formula is as follows: nz=2 h-1 Wherein nz represents the subdivision number of the underground depth under the maximum depth of the tree structure, h represents the maximum depth of the tree structure, and a binary tree structure under the maximum depth h is constructed;
building a tree structure initial model m when sampling times t=0 0 ＝(Τ k,j ,ρ 1 ,K,ρ k ) Where k is the number of nodes in the initial model, j is the depth level of the tree structure, j= (1, 2,., h), t k,j Representing the number of permutations of the tree structure of all different shapes composed of k node numbers at depth j, ρ k Representing conductivity parameters formed by the kth node;
building a markov chain starting from t=1, and sampling the existing tree structure m at the t-th time t In other words, the nodes of the tree structure can be classified into existing nodes, generatable nodes, and deceased nodes;
s3: optionally tree-like structure m t Tree structure change is carried out on the class of nodes of the model (a), and a new candidate model m' is formed; if the generatable node is selected, a new candidate tree structure node k '=k+1 is represented, and the new candidate model is m' = (t) k,j ,ρ 1 ,K,ρ k ,ρ k+1 ) The method comprises the steps of carrying out a first treatment on the surface of the If an erasable node is selected, a new candidate tree structure node k '=k+1 is represented, and the new candidate model is m' = (t) k,j ,ρ 1 ,K,ρ k-1 ) The method comprises the steps of carrying out a first treatment on the surface of the If one of the existing nodes is selected, it means that only the conductivity parameter ρ 'of the selected node is changed' n The new candidate model is m' = (t) k,j ,ρ 1 ,K,ρ′ n ,K,ρ k )；
S4: performing inverse transformation on the new candidate model m' under the tree structure obtained in the step S3 to obtain a suggested model of the conductivity structureAnd the advice model is +_ by the forward operator G>Mapping to a data space to obtain forward response;
s5: calculating a suggestion model based on the electromagnetic field response in S1 and the forward response obtained in S4Likelihood functions of (2);
s6: calculating based on the likelihood function in S5 to obtain the acceptance rate alpha of the suggested model, taking a random number between 0 and 1, judging whether the acceptance rate alpha is smaller than the random number, and accepting the new candidate model if yes, so as to obtain the existing tree structure when the t+1st sampling is as follows: m is m t+1 Otherwise reject the new candidate model, i.e. get the existing tree structure at the t+1st sampling as follows: m is m t+1 ＝m t ；
S7: repeating S3 to S6, and storing the obtained sampling model once every Z times of sampling until the sampling times t=Q times are finished;
s8: and (3) carrying out statistical calculation on all the sampling models obtained in the step (S7) to obtain the one-dimensional edge probability density distribution of the underground conductivity of the area.
Optionally, the inverse transformation in S4 specifically includes:
the scale function phi (a) is constructed to control the node currently existing, and the specific formula is as follows:
and for the kth node at tree depth j is defined as:
φ j,k (a)＝φ(2 j-1 a-e),e＝1,K,k；
wherein X is [0,1) (a) Expressed as a unit function over the interval 0, 1), any number a expressed as 0-1 has no meaning, j expressed as a level at which the node is located, and e expressed as a position of the level at which the node is located;
the decoupling function psi (a) is constructed to decouple the relation between the upper level node and the next level node, namely, the value of the father node of the range of the child node is covered by the value of the child node, and the specific formula is as follows:
similarly, for the kth node definition at tree depth j:
ψ j,k (a)＝ψ(2 j-1 a-e),e＝1,K,2 j ；
wherein Y is [0,1) (a) Represented as a unit function over interval 0, 1).
Optionally, the S5 includes:
s5.1: suggested model obtained according to inverse transformation of tree-structured candidate modelLeast square method is carried out on the forward calculation response of (2) and the electromagnetic field response obtained by measurement, and fitting difference +.>The specific formula is as follows:
wherein G represents the positive operator of the mapping of the model vector to the data space, C e Is a data covariance matrix, d is actually measured electromagnetic field response data;
s5.2: assuming that likelihood functions obey a standard normal distribution, likelihood functionsThe specific formula of (2) is as follows:
where n is the number of observations.
Optionally, the S6 includes:
s6.1: computing a priori distribution p (m t ) And a priori distribution p (m') of the new candidate model;
s6.2: calculating the existing tree structure model m in the t-th sampling t Inverse transformation of the resulting conductivity modelLikelihood function of +.>And a proposed model +.The new candidate tree-structure model m' at the time of the t-th sampling is inverse transformed>Likelihood function of +.>
S6.3: calculation from existing tree-structured model m t The proposed distribution probability q (m' |m) to the new candidate tree-structure model m t ) From the new candidate tree-structure model m' to the existing tree-structure model m t Is a proposed distribution q (m t |m′)；
S6.4: the acceptance rate alpha of the suggested model is calculated according to the prior distribution, the likelihood function and the suggested distribution calculated in the steps S6.1 to S6.3, and the specific formula is as follows:
wherein J is a jacobian matrix, |j|=1;
s6.5: taking a random number between 0 and 1, judging whether the acceptance rate alpha is smaller than the random number, and accepting a suggestion model, namely m if yes t+1 =m', otherwise reject the proposed model, i.e. m t+1 ＝m t 。
Optionally, the step S6.1 includes:
s6.1.1: calculating uniform distribution p (k) in a given node number range, wherein the specific formula is as follows:
wherein k is max Is the maximum node number k of the tree structure max ＝2 h -1；k min K is the minimum node number of the tree structure min Generally set to 1, indicating at least one node is owned;
s6.1.2: calculating the permutation number a priori p (T k I k), that is, the distribution of the number of permutation and combination of the tree structure in the case of a given node, the specific formula is as follows:
wherein T is k The number of permutations of the tree structure of all different shapes for the number of k nodes,is a standard binomial coefficient;
s6.1.3: calculating a priori distribution p (p) i |T k K), the specific formula is as follows:
s6.1.4: the prior distribution p (m) under the binary tree structure is calculated by adopting the same formula t ) And p (m'), wherein the a priori distribution p (m) t ) The specific formula of (2) is as follows:
optionally, the S6.3 includes:
s6.3.1: forming the current nodes on the current tree structure into a perturbable node group S v The generatable nodes form a generatable node group S b The resolvable nodes form a resolvable node group S d ；
S6.3.2: if the node selected in the S3 is the existing node, calculating disturbance suggestion distribution, if the node selected in the S3 is a generatable node, calculating generation suggestion distribution, and if the node selected in the S3 is an deceaseable node, calculating deceased suggestion distribution.
Optionally, the calculating of the disturbance recommendation distribution in S6.3.2 includes:
calculate a given S v Probability q (i|S) of selecting the ith node v ) The specific formula is as follows:
wherein |S v I is S v The number of elements;
calculating a probability q (Δρ) given the selection of conductivity values at the ith tree node i I), the specific formula is as follows:
wherein sigma is the standard deviation of the normal distribution of parameter disturbance, deltaρ i Selecting a conductivity change value on the ith tree node;
the specific formula of the disturbance recommendation distribution in S6.3.2 is as follows:
q(m′|m t )＝q(i|S v )q(Δρ i |i)；
the specific formula for generating the suggested distribution in S6.3.2 is as follows:
q(m′|m t )＝q(i|S b )q(ρ i |i)；
wherein q (ρ) i |i)＝p(ρ i |T k K), and q (i|S b ) The probability of selecting the location where the new node is placed is expressed as follows:
wherein, |S b The I is expressed as the number of the node sets which can be generated in the current model;
the calculation formula of the extinction proposal distribution in S6.3.2 is as follows:
wherein S' d Is a set of nodes that may be deleted after the proposed birth, |S '' d The i is expressed as the number of sets of resolvable nodes of the candidate model.
Optionally, in S6:
the expression of the proposed model acceptance rate α for the selected existing node is:
the expression of the proposed model acceptance rate α for the selection of generatable nodes is:
the expression of the proposed model acceptance rate α for the selection of the resolvable nodes is:
wherein, |S d The I is expressed as the number of the set of resolvable nodes of the current model, |S '' b The i is represented as the number of sets of nodes that can be generated for the candidate model.
Optionally, in S7, Z is more than or equal to 10 and less than or equal to 20, Q is not less than 1000000, and Z and Q are natural numbers.
In order to realize a physical space tree structure, the invention improves the tree structure of a wavelet domain, and in the original tree structure, the node coefficient value is expressed as the wavelet coefficient without any physical meaning, so that we cannot accurately give the prior distribution condition of related coefficients, which leads the conductivity value obtained by inversion to be possibly beyond the actual condition, therefore, the invention changes the wavelet domain tree structure into the physical space tree structure, and the node parameter of the tree structure is the prior distribution with the electrical distribution characteristic. According to the invention, a likelihood function of the suggestion model is obtained according to the actually measured underground electromagnetic field response and the forward modeling response calculation of the suggestion model, and the acceptance rate of the suggestion model is calculated, so that the suggestion model is screened to obtain a sampling model which is closer to the actual underground conductivity model, and the conductivity parameters of the sampling model can be changed by optional nodes, so that inversion parameters can be flexibly adjusted, and the sampling can be conveniently carried out for a plurality of times; the invention solves the problem that the prior distribution is difficult to select by the nodes of the wavelet domain tree structure by adopting the physical space tree structure, achieves the effects of automating parameters and multi-scale subdivision of inversion areas, saves the search of a Markov chain on a low probability area, effectively improves the sampling efficiency, can be extended to Gao Weibei leaf inversion, and provides an automatic parameterized thought for a Gao Weibei leaf inversion method.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings that are required in the embodiments or the description of the prior art will be briefly described, and it is obvious that the drawings in the following description are only some embodiments of the present invention, and other drawings may be obtained according to the structures shown in these drawings without inventive effort for a person skilled in the art.
FIG. 1 is a flow chart of a design method of a cross-dimensional Bayesian sampler based on a physical space tree structure in an embodiment of the present invention;
FIG. 2 is a schematic diagram of a sampling flow in an embodiment of the invention;
FIG. 3 is a diagram of three different types of nodes included in a tree structure used in the present invention;
FIG. 4 shows a different arrangement of the same nodes in the embodiment of the present invention, wherein: fig. 4 (a) to 4 (e) show all tree structure arrangements when tree nodes (black boxes are shown) are 3;
FIG. 5 is a diagram illustrating a one-dimensional electrical structure corresponding to a binary tree in accordance with an embodiment of the present invention, wherein: fig. 5 (a) shows the conductivity distribution of the underground uniform half space corresponding to the binary tree when only the root node exists, fig. 5 (b) shows the underground conductivity distribution corresponding to the case where the first two layers of nodes exist in the tree structure, and fig. 5 (c) shows the underground conductivity distribution corresponding to the case where the first three layers of nodes exist in the tree structure; z is Z min For a measured minimum subsurface depth, Z max For a measured maximum subsurface depth;
FIG. 6 is a diagram of mapping functions used in the conversion from a tree structure to a physical model in the practice of the present invention, wherein: fig. 6 (a), 6 (b), 6 (c) show graphs of scale functions at different scale, fig. 6 (d), 6 (e), 6 (f) show graphs of decoupling functions at different scale (the first number of subscripts indicates scale, the second number indicates scale);
FIG. 7 is a prior verification diagram of tree nodes in an embodiment of the present invention, wherein: fig. 7 (a), fig. 7 (b), and fig. 7 (c) are respectively represented as a tree structure in a binary tree, k max Results plots for λ=10, 20,30 for poisson coefficients=50;
FIG. 8 is a graph of a one-dimensional edge probability density distribution (the distribution of the true model is the black line) obtained by inversion in the practice of the present invention, wherein: the abscissa is the distribution of conductivity values, and the ordinate is the depth;
FIG. 9 is a graph of the distribution of RMS during sampling in an embodiment of the invention, wherein: the abscissa is the number of times of sampling, and the ordinate is the fitting difference of the data;
FIG. 10 is a graph of a sample node number distribution in an embodiment of the present invention, wherein: the abscissa is the node number, and the ordinate is the frequency;
fig. 11 is a result of sampling using a wavelet based tree structure, wherein: fig. 11 (a) shows 100 tens of thousands of times, fig. 11 (b) shows 500 tens of thousands of times, and fig. 11 (c) shows 1500 tens of thousands of times;
FIG. 12 is a graph showing the results of sampling using a physical space tree-based cross-dimensional Bayesian sampler in accordance with the present invention, in which: fig. 12 (a) shows the total number of samples as 100 ten thousand times, fig. 12 (b) shows the total number of samples as 500 ten thousand times, and fig. 12 (c) shows the result of the total number of samples as 1500 ten thousand times.
The achievement of the objects, functional features and advantages of the present invention will be further described with reference to the accompanying drawings, in conjunction with the embodiments.
Detailed Description
Embodiments of the invention are described in detail below with reference to the attached drawings, but the invention can be implemented in a number of different ways, which are defined and covered by the claims.
The invention provides a design method of a cross-dimensional Bayes sampler based on a physical space tree structure, and aims to solve the problems that prior distribution selection of nodes of a wavelet domain tree structure is difficult and sampling efficiency is low.
As shown in fig. 1, the design method of the cross-dimensional bayesian sampler based on the physical space tree structure comprises the following steps:
s1: measuring by using an electromagnetic measuring instrument to obtain the response of an underground electromagnetic field of a certain area, and determining the measured underground depth;
s2: constructing a tree structure initial model, and constructing a Markov chain, wherein data, an initial model m and prior information are input as shown in fig. 2; calculating initial model fitting difference, starting sampling, and from 1 to the maximum sampling times;
the details are as follows:
calculating the maximum depth of the tree structure, wherein the formula is as follows: nz=2 h-1 Wherein nz represents the subdivision number of the underground depth under the maximum depth of the tree structure, h represents the maximum depth of the tree structure, and a binary tree structure under the maximum depth h is constructed;
building a tree structure initial model m when sampling times t=0 0 ＝(Τ k,j ,ρ 1 ,K,ρ k ) Where k is the number of nodes in the initial model, j is the depth level of the tree structure, j= (1, 2,., h), t k,j Representing the number of permutations of the tree structure of all different shapes composed of k node numbers at depth j, ρ k The conductivity parameter formed by the kth node is shown in fig. 5;
building a markov chain starting from t=1, and sampling the existing tree structure m at the t-th time t In other words, the nodes of the tree structure can be classified into existing nodes, generatable nodes, and deceased nodes;
s3: optionally tree-like structure m t Tree structure change is carried out on the class of nodes of the model (a), and a new candidate model m' is formed; if the generatable node is selected, a new candidate tree structure node k '=k+1 is represented, and the new candidate model is m' = (t) k,j ,ρ 1 ,K,ρ k ,ρ k+1 ) The method comprises the steps of carrying out a first treatment on the surface of the If an erasable node is selected, a new candidate tree structure node k '=k+1 is represented, and the new candidate model is m' = (t) k,j ,ρ 1 ,K,ρ k-1 ) The method comprises the steps of carrying out a first treatment on the surface of the If one of the existing nodes is selected, it means that only the conductivity parameter ρ 'of the selected node is changed' n The new candidate model is m' = (t) k,j ,ρ 1 ,K,ρ′ n ,K,ρ k )；
S4: performing inverse transformation on the new candidate model m' under the tree structure obtained in the step S3 to obtain a suggested model of the conductivity structureAnd will be calculated by the positive operator GAdvice model->Mapping to a data space to obtain forward response;
specifically, as shown in fig. 6, the inverse transformation in S4 specifically includes:
the scale function phi (a) is constructed to control the node currently existing, and the specific formula is as follows:
and for the kth node at tree depth j, it is possible to define:
φ j,k (a)＝φ(2 j-1 a-e),e＝1,K,k；
wherein X is [0,1) (a) Expressed as a unit function over the interval [0,1 ], a expressed as any number from 0 to 1 has no specific meaning, j expressed as a level at which the node is located, and e expressed as a position of the level at which the node is located;
the decoupling function psi (a) is constructed to decouple the relation between the upper level node and the next level node, namely, the value of the father node of the range of the child node is covered by the value of the child node, and the specific formula is as follows:
similarly, for the kth node definition at tree depth j:
ψ j,k (a)＝ψ(2 j-1 a-e),e＝1,K,2 j ；
wherein Y is [0,1) (a) Represented as a unit function over interval 0, 1).
S5: calculating a suggestion model based on the electromagnetic field response in S1 and the forward response obtained in S4Likelihood functions of (2);
the method specifically comprises the following steps:
s5.1: suggested model obtained according to inverse transformation of tree-structured candidate modelLeast square method is carried out on the forward response of the model (A) and the electromagnetic field response obtained by measurement, and fitting difference +.>The specific formula is as follows:
wherein G represents the positive operator of the mapping of the model vector to the data space, C e Is a data covariance matrix, d is actually measured electromagnetic field response data;
s5.2: assuming that likelihood functions obey a standard normal distribution, likelihood functionsThe distribution of (2) is specifically expressed as follows:
where n is the number of observations.
S6: calculating based on the likelihood function in S5 to obtain the acceptance rate alpha of the suggested model, taking a random number between 0 and 1, judging whether the acceptance rate alpha is smaller than the random number, and accepting the new candidate model if yes, so as to obtain the existing tree structure when the t+1st sampling is as follows: m is m t+1 Otherwise reject the new candidate model, i.e. get the existing tree structure at the t+1st sampling as follows: m is m t+1 ＝m t ；
The method specifically comprises the following steps:
s6.1: computing a priori distribution p (m t ) And a priori distribution p (m') of the new candidate model;
s6.2: calculating the existing tree structure model m in the t-th sampling t Inverse transformation of the resulting conductivity modelLikelihood function of +.>And a proposed model +.The new candidate tree-structure model m' at the time of the t-th sampling is inverse transformed>Likelihood function of +.>
S6.3: calculation from existing tree-structured model m t The proposed distribution probability q (m' |m) to the new candidate tree-structure model m t ) From the new candidate tree-structure model m' to the existing tree-structure model m t Is a proposed distribution q (m t |m′)；
S6.4: the acceptance rate alpha of the suggested model is calculated according to the prior distribution, the likelihood function and the suggested distribution calculated in the steps S6.1 to S6.3, and the specific formula is as follows:
wherein J is a jacobian matrix, |j|=1;
s6.5: taking a random number between 0 and 1, judging whether the acceptance rate alpha is smaller than the random number, and accepting a suggestion model, namely m if yes t+1 =m', otherwise reject the proposed model, i.e. m t+1 ＝m t 。
S6.1 includes:
s6.1.1: calculating uniform distribution p (k) in a given node number range, wherein the specific formula is as follows:
wherein k is max Is the maximum node number k of the tree structure max ＝2 h -1；k min K is the minimum node number of the tree structure min Generally set to 1 means at least one node is owned;
s6.1.2: calculating the permutation number a priori p (T k I k), that is, the distribution of the number of permutation and combination of the tree structure in the case of a given node, the specific formula is as follows:
wherein T is k The number of permutations of the tree structure of all different shapes for the number of k nodes,is a standard binomial coefficient;
in one embodiment, given the number of nodes k=3, the tree structure may be arranged as shown in fig. 4. S6.1.3: calculating a priori distribution p (p) i |T k K), the specific formula is as follows:
s6.1.4: the prior distribution p (m) under the binary tree structure is calculated by adopting the same formula t ) And p (m'), wherein the a priori distribution p (m) t ) The specific formula of (2) is as follows:
s6.3 comprises:
s6.3.1: forming the current nodes on the current tree structure into a perturbable node group S v The generatable nodes form a generatable node groupS b The resolvable nodes form a resolvable node group S d ；
As shown in FIG. 3, the resolvable node group represents all nodes of the child nodes which are not active on the current tree structure, the perturbable node group represents all nodes of the child nodes which are not active on the current tree structure and only one child node, and the generatable node group is a collection of hollow nodes in the tree structure.
S6.3.2: if the node selected in the S3 is the existing node, calculating disturbance suggestion distribution, if the node selected in the S3 is a generatable node, calculating generation suggestion distribution, and if the node selected in the S3 is an deceaseable node, calculating deceased suggestion distribution.
The calculation of the disturbance recommendation distribution in S6.3.2 includes:
calculate a given S v Probability q (i|S) of selecting the ith node v ) The specific formula is as follows:
wherein |S v I is S v The number of elements;
calculating a probability q (Δρ) given the selection of conductivity values at the ith tree node i I), the specific formula is as follows:
wherein sigma is the standard deviation of the normal distribution of parameter disturbance, deltaρ i Selecting a conductivity change value on the ith tree node;
the specific formula for calculating the disturbance recommendation distribution in S6.3.2 is as follows:
q(m′|m t )＝q(i|S v )q(Δρ i |i)；
the specific formula for generating the suggested distribution in S6.3.2 is as follows:
q(m′|m t )＝q(i|S b )q(ρ i |i)；
wherein q (ρ) i |i)＝p(ρ i |T k K), and q (i|S b ) The probability of selecting the location where the new node is placed is expressed as follows:
wherein, |S b The I is expressed as the number of the node sets which can be generated in the current model;
the calculation formula of the extinction proposal distribution in S6.3.2 is as follows:
wherein S' d Is a set of nodes that may be deleted after the proposed birth, |S '' d The i is expressed as the number of sets of resolvable nodes of the candidate model.
S5:
the expression of the proposed model acceptance rate α for the selected existing node is:
the expression of the proposed model acceptance rate α for the selection of generatable nodes is:
the expression of the proposed model acceptance rate α for the selection of the resolvable nodes is:
wherein, |S d The I is expressed as the number of the set of resolvable nodes of the current model, |S '' b The i is represented as the number of sets of nodes that can be generated for the candidate model.
S7: let t=t+1, repeat S3 to S6, and save the obtained sampling model once every Z samplings until the sampling times t=q ends. Wherein Z is more than or equal to 10 and less than or equal to 20, Q is not less than 1000000, and Z and Q are natural numbers;
s8: and (3) carrying out statistical calculation on all the sampling models obtained in the step (S7), and obtaining the one-dimensional edge probability density distribution of the underground conductivity of the area, as shown in figure 8.
For such a tree-structured cross-dimensional sampler, it is critical whether the M-H acceptance criteria constructed above are correct. While a key check on the correctness of the acceptance criteria for a tree-structured cross-dimensional sampler is to test whether there is a posterior deviation of the criterion for the number of nodes k. One effective method is that the sampler operates a large number of steps, and regards likelihood function items as constant items, and finally compares the distribution condition of the sampled node number with the estimated distribution condition to judge whether the distribution accords with the expectations. For the likelihood function, a constant is maintained, so the perturbation step is not performed in the cross-dimensional sampler, the likelihood ratio in the generation step and the extinction step is kept constant at 1, and the probability of performing p (birth) =p (delay) =0.5.
In this embodiment, the number of nodes k ε [1, k ] max ]100 ten thousand samples were tested, k max For a manually-set maximum node value, k max =50. Since λ is in a poisson distribution, a priori p about the number of nodes k can be obtained before sampling prior (k) The following is shown:
the result of setting λ to 10,20, and 30 is shown in fig. 7. The black bar graph in the graph represents the distribution of the number of nodes after sampling, and the black line represents the prior distribution to which the nodes should obey. From the figure it can be seen that the sampled distributions at λ=10, λ=20 and λ=30 are quite identical to the a priori known distributions.
In order to verify inversion effect of the method provided by the invention, thus proving that the required uncertainty information can be obtained, performing model inspection, designing an underground conductivity model shown in table 1, sampling by a designed sampler to obtain a posterior probability model, and describing in detail as follows:
table 1: conductivity and interface depth information
Layer(s) 1 2 3 4 5
σ(S/m) 0.001 0.01 0.1 0.0025 0.01
z(m) 1000 1000 1000 1000 ∞
The cross-dimension sampler designed by the invention is realized by utilizing Fortran language programming, and a personal notebook computer used for running a program is configured as follows: CPU-InterCore i7-8700, main frequency 3.2GHz, running memory 8.00GB. FIG. 8 is a graph of a one-dimensional edge probability density distribution of a distribution of an underground medium sampled by a sampler according to the present invention, wherein the shade of the edge probability density distribution indicates the conductivity probability of a certain region at that depth, and the darker the shade indicates the higher the probability. From the results obtained by inversion, the 5-layer model is better recovered. The first, second, third and last layers have better edge probability density constraint and small uncertainty, so that the sampling efficiency of the sampler for the layers is high, the uncertainty range of the fourth layer is large, the fact that the magnetotelluric method is insensitive to high resistance and the expressive force of the magnetotelluric method is supposed to have physical characteristics is shown, but from the probability of the magnetotelluric method, the high probability area is better corresponding to the real situation of the model, and the method can well extract the characteristic response of the high resistance layer in the data. And the sampling by the samplers shown in fig. 9 and 10 has good mixing and constraint on the number of parameters.
When inverting a conductivity model of three subsurface layers, when sampling using a wavelet based tree structure, it can be seen that the resulting distribution in FIG. 11 does not match the true model distribution until 500 tens of thousands of times after sampling. When the cross-dimensional Bayesian sampler based on the physical space tree structure is used for sampling, the distribution of the real model can be searched by only carrying out one million times of sampling in the figure 12. Therefore, the sampling efficiency is remarkably improved.
The foregoing description is only of the preferred embodiments of the present invention and is not intended to limit the scope of the invention, and all equivalent structural changes made by the description of the present invention and the accompanying drawings or direct/indirect application in other related technical fields are included in the scope of the invention.

Claims (9)

1.一种基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，该设计方法包括：1. A design method for a cross-dimensional Bayesian sampler based on the tree structure of physical space, which is characterized in that the design method includes: S1：采用电磁测量仪器测量获得某区域的地下电磁场响应，并确定测量的地下深度；S1: Use electromagnetic measuring instruments to measure and obtain the underground electromagnetic field response in a certain area, and determine the measured underground depth; S2：构建树状结构初始模型，构建马尔科夫链，具体是：S2: Construct an initial model of tree structure and build a Markov chain, specifically: 计算树状结构的最大深度，公式为：nz＝2h-1，其中，h表示树状结构的最大深度，nz表示树状结构最大深度下的地下深度的剖分数量，并构建最大深度h下的二叉树结构；Calculate the maximum depth of the tree structure, the formula is: nz=2 h-1 , where h represents the maximum depth of the tree structure, nz represents the number of subdivisions of the underground depth under the maximum depth of the tree structure, and constructs the maximum depth h The binary tree structure below; 构建采样次数t＝0时的树状结构初始模型m0＝(Τk,j,ρ1,K,ρk)，其中，k为初始模型中的节点数，j为树状结构的深度等级，j＝(1,2,...,h)，Τk,j表示为深度j下的k个节点数所构成全部不同形状的树状结构的排列数，ρk表示第k个节点所构成的电导率参数；Construct the initial model of the tree structure m 0 = (Τ k,j ,ρ 1 ,K,ρ k ) when the sampling number t=0, where k is the number of nodes in the initial model and j is the depth level of the tree structure , j=(1,2,...,h), Τ k,j represents the number of arrangements of all tree-like structures of different shapes composed of k nodes at depth j, ρ k represents the k node The constituted conductivity parameter; 从t＝1时开始构建马尔科夫链，对第t次采样时的现有树状结构mt而言，将树状结构的节点分类为现有节点、可生成节点以及可消亡节点；A Markov chain is constructed starting from t=1. For the existing tree structure m t at the t-th sampling time, the nodes of the tree structure are classified into existing nodes, generateable nodes and perishable nodes; S3：任选树状结构mt的一类节点进行树状结构变化，形成新候选模型m′；若选择可生成节点，则表示新候选树状结构节点k'＝k+1，新候选模型为m′＝(Τk,j,ρ1,K,ρk,ρk+1)；若选择可消亡节点，则表示新候选树状结构节点k'＝k-1，新候选模型为m′＝(Τk,j,ρ1,K,ρk-1)；若选择现有节点中的一个时表示仅仅改变选中节点的电导率参数ρ′n，新候选模型为m′＝(Τk,j,ρ1,K,ρ′n,K,ρk)；S3: Select a type of node in the tree structure m t to change the tree structure to form a new candidate model m′; if a node can be generated, it means that the new candidate tree structure node k′=k+1, the new candidate model is m′=(Τ k,j ,ρ 1 ,K,ρ k ,ρ k+1 ); if a perishable node is selected, it means that the new candidate tree structure node k'=k-1, and the new candidate model is m ′=(Τ k,j ,ρ 1 ,K,ρ k-1 ); if selecting one of the existing nodes means only changing the conductivity parameter ρ′ n of the selected node, the new candidate model is m′=(Τ k,j ,ρ 1 ,K,ρ′ n ,K,ρ k ); S4：将S3中得到的树状结构下的新候选模型m′进行逆变换得到电导率结构的建议模型并通过正演算子G将建议模型/>映射到数据空间得到正演响应；S4: Inversely transform the new candidate model m′ under the tree structure obtained in S3 to obtain the recommended model of conductivity structure And pass the forward operator G to propose the model/> Map to data space to obtain forward response; S5：基于S1中的电磁场响应和S4中得到的正演响应计算建议模型的似然函数；S5: Calculate the proposed model based on the electromagnetic field response in S1 and the forward response obtained in S4 the likelihood function; S6：基于S5中的似然函数计算得到建议模型的接受率α，判断接受率α是否小于此随机数，若是则接受新候选模型，即得到第t+1次采样时的现有树状结构为：mt+1＝m′，否则拒绝新候选模型，即得到第t+1次采样时的现有树状结构为：mt+1＝mt；S6: Calculate the acceptance rate α of the proposed model based on the likelihood function in S5, and determine whether the acceptance rate α is less than this random number. If so, accept the new candidate model, that is, obtain the existing tree structure at the t+1th sampling. is: m t+1 =m′, otherwise the new candidate model is rejected, that is, the existing tree structure at the t+1th sampling is obtained: m t+1 =m t ; S7：令t＝t+1，重复S3至S6，并每进行Z次采样保存一次得到的采样模型，直至采样次数t＝Q次结束；S7: Let t=t+1, repeat S3 to S6, and save the obtained sampling model every Z sampling, until the number of sampling t=Q times; S8：将S7中保存得到的所有采样模型进行统计计算，得到此区域的地下电导率的一维边缘概率密度分布。S8: Perform statistical calculations on all sampling models saved in S7 to obtain the one-dimensional marginal probability density distribution of underground conductivity in this area. 2.如权利要求1所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S4中的逆变换具体包括：2. The design method of a cross-dimensional Bayesian sampler based on the tree structure of physical space as claimed in claim 1, characterized in that the inverse transformation in S4 specifically includes: 构建尺度函数φ(a)，以控制当前存在的节点，具体公式如下：Construct a scale function φ(a) to control the currently existing nodes. The specific formula is as follows: 而对于树状结构深度为j下的第k个节点定义为：The k-th node under the tree structure depth j is defined as: φj,k(a)＝φ(2j-1a-e),e＝1,K,k；φ j,k (a)=φ(2 j-1 ae),e=1,K,k; 其中：X[0,1)(a)表示为在区间[0,1)上的单位函数；a表示为0-1上的任意数无具体意义；j表示为节点所在的等级；e表示为节点所在该等级的位置；Among them : The location of the node at this level; 构建解耦函数ψ(a)，以解耦上一等级节点与其下一等级节点的关系，即用其子节点的值覆盖子节点所在范围的父节点的值，具体公式如下：Construct a decoupling function ψ(a) to decouple the relationship between the upper level node and the next level node, that is, use the value of its child node to cover the value of the parent node in the range of the child node. The specific formula is as follows: 同理，对于树状结构深度为j下的第k个节点定义为：In the same way, the k-th node with a tree structure depth of j is defined as: ψj,k(a)＝ψ(2j-1a-e),e＝1,K,2j；ψ j,k (a)=ψ(2 j-1 ae),e=1,K,2 j ; 其中，Y[0,1)(a)表示为在区间[0,1)上的单位函数。Among them, Y [0,1) (a) is expressed as a unit function on the interval [0,1). 3.如权利要求1所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S5包括：3. The design method of a cross-dimensional Bayesian sampler based on the tree structure of physical space as claimed in claim 1, characterized in that the S5 includes: S5.1：根据树状结构候选模型逆变换得到的建议模型的正演响应和测量得到的电磁场响应进行最小二乘法计算，得到其拟合差/>具体公式如下：S5.1: Recommended model obtained by inverse transformation of candidate models based on tree structure The forward modeling response and the measured electromagnetic field response are calculated by the least square method to obtain the fitting difference/> The specific formula is as follows: 其中，G代表的是模型向量映射到数据空间中的正演算子，Ce是数据协方差矩阵，d为实际测量到的电磁场响应数据；Among them, G represents the forward operator that maps the model vector to the data space, C e is the data covariance matrix, and d is the actual measured electromagnetic field response data; S5.2：假设似然函数服从标准正态分布，似然函数的具体公式如下：S5.2: Assuming that the likelihood function obeys the standard normal distribution, the likelihood function The specific formula is as follows: 其中，n为观测数据的数量。Among them, n is the number of observation data. 4.如权利要求3所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S6包括：4. The design method of a cross-dimensional Bayesian sampler based on the tree structure of physical space as claimed in claim 3, characterized in that the S6 includes: S6.1：计算给定二叉树结构下的当前现有树状结构的先验分布p(mt)以及新候选模型的先验分布p(m′)；S6.1: Calculate the prior distribution p(m t ) of the current existing tree structure under the given binary tree structure and the prior distribution p(m′) of the new candidate model; S6.2：计算第t次采样时的现有树状结构模型mt逆变换得到的电导率模型的似然函数以及第t次采样时的新候选树状结构模型m′经过逆变换得到的建议模型/>的似然函数/> S6.2: Calculate the conductivity model obtained by the inverse transformation of the existing tree structure model m t at the t-th sampling the likelihood function of And the suggested model obtained by inverse transformation of the new candidate tree structure model m′ at the t-th sampling/> Likelihood function/> S6.3：计算从现有树状结构模型mt到新候选树状结构模型m′的建议分布概率q(m′|mt)以及从新候选树状结构模型m′到现有树状结构模型mt的建议分布q(mt|m′)；S6.3: Calculate the proposed distribution probability q(m′|m t ) from the existing tree structure model m t to the new candidate tree structure model m′ and from the new candidate tree structure model m′ to the existing tree structure Proposed distribution q(m t |m′) of model m t ; S6.4：根据S6.1至S6.3中计算得到的先验分布、似然函数和建议分布计算建议模型的接受率α，具体公式如下：S6.4: Calculate the acceptance rate α of the proposed model based on the prior distribution, likelihood function and suggested distribution calculated in S6.1 to S6.3. The specific formula is as follows: 其中，J为雅克比矩阵，|J|＝1；Among them, J is the Jacobian matrix, |J|=1; S6.5：取0-1之间的随机数，判断接受率α是否小于此随机数，若是则接受建议模型，即mt+1＝m′，否则拒绝建议模型，即mt+1＝mt。S6.5: Take a random number between 0-1 and determine whether the acceptance rate α is less than this random number. If so, accept the suggested model, that is, m t+1 = m′, otherwise reject the suggested model, that is, m t+1 = mt . 5.如权利要求4所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S6.1中包括：5. The design method of a cross-dimensional Bayesian sampler based on the tree structure of physical space as claimed in claim 4, characterized in that the S6.1 includes: S6.1.1：计算给定节点数范围下的均匀分布p(k)，具体公式如下：S6.1.1: Calculate the uniform distribution p(k) under a given range of node numbers. The specific formula is as follows: 其中，kmax为树状结构的最大节点数，kmax＝2h-1；kmin为树状结构的最小节点数，kmin一般设置为1，表示至少拥有一个节点；Among them, k max is the maximum number of nodes in the tree structure, k max =2 h -1; k min is the minimum number of nodes in the tree structure, and k min is generally set to 1, indicating that there is at least one node; S6.1.2：计算给定节点数下的排列数先验p(Tk|k)，即为在给定节点的情况下，树状结构的排列组合数的分布情况，具体公式如下：S6.1.2: Calculate the number of permutations prior p(T k |k) for a given number of nodes, which is the distribution of the number of permutations and combinations of the tree structure given the number of nodes. The specific formula is as follows: 其中，Tk为k个节点数所构成全部不同形状的树状结构的排列数， 为标准二项式系数；Among them, T k is the number of arrangements of all tree structures of different shapes composed of k nodes, is the standard binomial coefficient; S6.1.3：计算给定每个节点中参数的先验分布p(ρi|Tk,k)，具体公式如下：S6.1.3: Calculate the prior distribution p(ρ i |T k ,k) given the parameters in each node. The specific formula is as follows: S6.1.4：采用同一公式分别计算二叉树结构下的先验分布p(mt)和p(m′)，其中先验分布p(mt)的具体公式如下：S6.1.4: Use the same formula to calculate the prior distribution p(m t ) and p(m′) under the binary tree structure. The specific formula of the prior distribution p(m t ) is as follows: 6.如权利要求5所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S6.3包括：6. The method for designing a cross-dimensional Bayesian sampler based on a physical space tree structure according to claim 5, wherein S6.3 comprises: S6.3.1：将当前树状结构上的现有节点构成可扰动节点组Sv，可生成节点构成可生成节点组Sb，可消亡节点构成可消亡节点组Sd；S6.3.1: Configure the existing nodes on the current tree structure into a perturbed node group S v , the generateable nodes into a generateable node group S b , and the perishable nodes into a perishable node group S d ; S6.3.2：若S3中选择的节点为现有节点，则计算扰动建议分布，若S3中选择的节点为可生成节点，则计算生成建议分布，若S3中选择的节点为可消亡节点，则计算消亡建议分布。S6.3.2: If the node selected in S3 is an existing node, calculate the perturbation suggestion distribution. If the node selected in S3 is a generative node, then calculate the recommended distribution. If the node selected in S3 is a perishable node, then Calculate the distribution of death suggestions. 7.如权利要求6所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S6.3.2中扰动建议分布的计算包括：7. The design method of a cross-dimensional Bayesian sampler based on a physical space tree structure as claimed in claim 6, characterized in that the calculation of the perturbation suggestion distribution in S6.3.2 includes: 计算给定Sv选择第i个节点的概率q(i|Sv)，具体公式如下：Calculate the probability q(i|S v ) of selecting the i-th node given S v . The specific formula is as follows: 其中|Sv|是Sv的元素个数；Where |S v | is the number of elements of S v ; 计算给定在第i个树节点上选择电导率值概率q(Δρi|i)，具体公式如下：Calculate the probability q(Δρ i |i) of selecting the conductivity value given on the i-th tree node. The specific formula is as follows: 其中，σ为参数扰动的正态分布的标准差，Δρi第i个树节点上的电导率变化值；Among them, σ is the standard deviation of the normal distribution of parameter perturbation, and Δρ is the conductivity change value on the i-th tree node; 所述S6.3.2中扰动建议分布的具体公式如下：The specific formula for the proposed distribution of disturbances in S6.3.2 is as follows: q(m′|mt)＝q(i|Sv)q(Δρi|i)；q(m′|m t )=q(i|S v )q(Δρ i |i); 所述S6.3.2中生成建议分布的具体公式如下：The specific formula for generating the suggested distribution in S6.3.2 is as follows: q(m′|mt)＝q(i|Sb)q(ρi|i)；q(m′|m t )=q(i|S b )q(ρ i |i); 其中，q(ρi|i)＝p(ρi|Tk,k),且q(i|Sb)表示选择放置新节点的位置的概率，具体公式如下：Among them, q(ρ i |i)=p(ρ i |T k ,k), and q(i|S b ) represents the probability of choosing the location to place the new node. The specific formula is as follows: 其中，|Sb|表示为当前模型的可生成节点集合数量；Among them, |S b | represents the number of node sets that can be generated in the current model; 所述S6.3.2中消亡建议分布的计算公式为：The calculation formula for the distribution of extinction suggestions in S6.3.2 is: 其中，S′d是在提议的出生后可能被删除的节点集，|S′d|表示为候选模型的可消亡节点集合数量。Among them, S′ d is the set of nodes that may be deleted after the birth of the proposal, and |S′ d | represents the number of extinct nodes in the candidate model. 8.如权利要求7所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S6中：8. The design method of a cross-dimensional Bayesian sampler based on the tree structure of physical space as claimed in claim 7, characterized in that, in the S6: 选择现有节点的建议模型接受率α的表达式为：The expression for selecting the proposed model acceptance rate α of existing nodes is: 选择可生成节点的建议模型接受率α的表达式为：The expression of the proposed model acceptance rate α for selecting the node that can be generated is: 选择可消亡节点的建议模型接受率α的表达式为：The expression of the recommended model acceptance rate α for selecting extinct nodes is: 其中，|Sd|表示为当前模型的可消亡节点集合数量，|Sb′|表示为候选模型的可生成节点集合数量。Among them, |S d | represents the number of perishable node sets of the current model, and |S b ′| represents the number of generative node sets of the candidate model. 9.如权利要求1至8中任一项所述的基于物理空间树状结构的跨维贝叶斯采样器的设计方法，其特征在于，所述S7中，10≤Z≤20，Q不小于1000000，Z、Q取自然数。9. The design method of a cross-dimensional Bayesian sampler based on a tree structure of physical space according to any one of claims 1 to 8, characterized in that in the S7, 10≤Z≤20, Q does not If it is less than 1,000,000, Z and Q are natural numbers.
CN202410078165.4A 2024-01-19 2024-01-19 Design method of cross-dimension Bayesian sampler based on physical space tree structure Active CN117574790B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202410078165.4A CN117574790B (en) 2024-01-19 2024-01-19 Design method of cross-dimension Bayesian sampler based on physical space tree structure

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202410078165.4A CN117574790B (en) 2024-01-19 2024-01-19 Design method of cross-dimension Bayesian sampler based on physical space tree structure

Publications (2)

Publication Number Publication Date
CN117574790A CN117574790A (en) 2024-02-20
CN117574790B true CN117574790B (en) 2024-03-26

Family

ID=89864865

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202410078165.4A Active CN117574790B (en) 2024-01-19 2024-01-19 Design method of cross-dimension Bayesian sampler based on physical space tree structure

Country Status (1)

Country Link
CN (1) CN117574790B (en)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN119986830B (en) * 2025-04-15 2025-06-10 中南大学 Magnetotelluric probability inversion method based on image compression

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2010151354A1 (en) * 2009-06-26 2010-12-29 Exxonmobil Upstream Research Company Constructing resistivity models from stochastic inversion
CN113177330A (en) * 2021-05-27 2021-07-27 吉林大学 Transient electromagnetic rapid statistical inversion method
CN113553773A (en) * 2021-08-16 2021-10-26 吉林大学 Inversion method of ground-air electromagnetic data based on Bayesian framework combined with neural network
CN114065511A (en) * 2021-11-15 2022-02-18 中南大学 Numerical simulation method, device, equipment and medium for two-dimensional forward modeling of magnetotelluric under undulating terrain
CN114896564A (en) * 2022-05-23 2022-08-12 武汉市市政建设集团有限公司 Transient electromagnetic two-dimensional Bayesian inversion method adopting self-adaptive Thiessen polygon parameterization
CN116992668A (en) * 2023-08-03 2023-11-03 中国科学院地质与地球物理研究所 Transient electromagnetic self-adaptive transverse constraint inversion method

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110869814B (en) * 2017-07-06 2022-06-14 雪佛龙美国公司 System and method for full waveform inversion of seismic data

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2010151354A1 (en) * 2009-06-26 2010-12-29 Exxonmobil Upstream Research Company Constructing resistivity models from stochastic inversion
CN113177330A (en) * 2021-05-27 2021-07-27 吉林大学 Transient electromagnetic rapid statistical inversion method
CN113553773A (en) * 2021-08-16 2021-10-26 吉林大学 Inversion method of ground-air electromagnetic data based on Bayesian framework combined with neural network
CN114065511A (en) * 2021-11-15 2022-02-18 中南大学 Numerical simulation method, device, equipment and medium for two-dimensional forward modeling of magnetotelluric under undulating terrain
CN114896564A (en) * 2022-05-23 2022-08-12 武汉市市政建设集团有限公司 Transient electromagnetic two-dimensional Bayesian inversion method adopting self-adaptive Thiessen polygon parameterization
CN116992668A (en) * 2023-08-03 2023-11-03 中国科学院地质与地球物理研究所 Transient electromagnetic self-adaptive transverse constraint inversion method

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于小波树状结构的二维大地电磁跨维贝叶斯反演;田圣琦;中国优秀硕士学位论文全文数据库 基础科学辑;20240115(第01期);A011-301 *

Also Published As

Publication number Publication date
CN117574790A (en) 2024-02-20

Similar Documents

Publication Publication Date Title
Mamalakis et al. Neural network attribution methods for problems in geoscience: A novel synthetic benchmark dataset
Maschio et al. Probabilistic history matching using discrete Latin Hypercube sampling and nonparametric density estimation
CN114154427B (en) Volume fracturing fracture expansion prediction method and system based on deep learning
Jiang et al. Electrical resistivity imaging inversion: An ISFLA trained kernel principal component wavelet neural network approach
Cherkassky et al. Computational intelligence in earth sciences and environmental applications: Issues and challenges
Ma et al. Integration of artificial intelligence and production data analysis for shale heterogeneity characterization in steam-assisted gravity-drainage reservoirs
CN103336305B (en) A kind of method dividing Sandstone Gas Reservoir high water cut based on gray theory
CN114742325B (en) A method and system for predicting ground settlement during step-method construction of subway tunnels
CN117574790B (en) Design method of cross-dimension Bayesian sampler based on physical space tree structure
Loquin et al. A fuzzy interval analysis approach to kriging with ill-known variogram and data
CN105607122A (en) Seismic texture extraction and enhancement method based on total variation seismic data decomposition model
CN119888104A (en) Method, medium and system for building geometric model of overlying strata structure
CN119293550A (en) A classification and prediction method based on InSAR time series large-scale deformation results
CN117787054A (en) Earthquake travel time and magnetotelluric sounding joint inversion method
Qiu et al. Novel multi-spatial receptive field (MSRF) XGBoost method for predicting geological cross-section based on sparse borehole data
Moghaddam et al. Application of machine learning methods in inferring surface water groundwater exchanges using high temporal resolution temperature measurements
CN114881171B (en) Land shale oil lithology type identification method and system based on convolutional neural network
Guo et al. Groundwater depth forecasting using configurational entropy spectral analyses with the optimal input
Garini et al. Filling-well: An effective technique to handle incomplete well-log data for lithology classification using machine learning algorithms
Abdulkhaleq et al. Advanced machine learning for missing petrophysical property imputation applied to improve the characterization of carbonate reservoirs
Shao et al. Flexible modeling of nonstationary extremal dependence using spatially fused LASSO and ridge penalties
Hu et al. Investigate the rainfall-runoff relationship and hydrological concepts inside LSTM
Gallop Facies probability from mixture distributions with non-stationary impedance errors
CN116224456A (en) Transient electromagnetic inversion method based on mixed density network
CN117388933A (en) Neural network nuclear magnetic logging curve prediction method and device based on feature enhancement

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