收藏切换
GPR-based model validation method for small samples
收藏切换
PDF
Fan YANG1, 2, Ping MA1, 2, Huan ZHANG1, 2, Wei LI1, 2, *, Ming YANG1, 2
Journal of Systems Engineering and Electronics | 2026, 37(3) : 867 - 877
Less
收藏切换
Journal of Systems Engineering and Electronics | 2026, 37(3): 867-877
SYSTEMS ENGINEERING
GPR-based model validation method for small samples
Full
Fan YANG1, 2, Ping MA1, 2, Huan ZHANG1, 2, Wei LI1, 2, *, Ming YANG1, 2
Affiliations
  • 1Control and Simulation Center, Harbin Institute of Technology, Harbin 150080, China
  • 2National Key Laboratory of Complex System Modeling and Simulation, Harbin 150080, China
Published: 2026-06-18 doi: 10.23919/JSEE.2026.000113
Outline
收藏切换

Validation for simulation models often confronts challenges with small samples due to the costs of time and money. To address this issue, this paper presents a validation method for small-sample dynamic outputs based on Gaussian process regression (GPR) models. Firstly, a validation framework based on Bayes statistics is proposed, shifting the focus from merely analyzing validation data to a more comprehensive analysis of posterior distributions. Subsequently, the posterior distributions of both the simulation outputs and the reference data are separately captured through segmented GPR. Then, the consistency of these posterior distributions is evaluated in terms of the central tendency and the distribution range. This consistency serves as a quantitative measure of the simulation model’s credibility, expressed as a value ranging from 0 to 1, where a value closer to 1 indicates higher credibility. Finally, the effectiveness of this validation method is demonstrated through a numerical example and an application example, highlighting its capability in uncertainty description and adaptability to small samples.

model validation  /  small samples  /  Gaussian process regression (GPR)  /  Bayes statistics
Fan YANG, Ping MA, Huan ZHANG, Wei LI, Ming YANG. GPR-based model validation method for small samples[J]. Journal of Systems Engineering and Electronics, 2026 , 37 (3) : 867 -877 . DOI: 10.23919/JSEE.2026.000113
Nowadays, simulation technology has entered into a new stage focused on complex system simulations [1]. Due to the complexity of the simulation models, issues such as lack of reference data, coupling between sub-models, uncertainty and variety of outputs have raised concerns about model validation [2]. Validating simulation results is a common method for assessing the credibility of the simulation models by measuring the consistency between the simulation outputs and the reference data [3,4].
To validate the simulation results under uncertainty, it is necessary to analyze the consistency of two sets of multi-sample validation data. Various methods have been proposed to solve this problem, including validation methods based on copulas [5], Bayesian networks [6], differences in multivariate probability distributions [7], grey cloud clustering [8], and improved Bayesian factor [9]. However, due to prolonged experimental cycles, substantial costs, and intricate processes, the outputs of both real systems and simulation models are typically small samples, which makes these methods unsuitable. For instance, in experiments such as airplane flights and rocket launches, the sample size of validation data is significantly limited.
Small-sample data is defined as the data with a sample size not exceeding 30 [10]. Various scholars have explored simulation validation techniques designed for small samples, but the broader practical application of these methods remains somewhat limited. The validation method based on bootstrap and Bayes parameter estimation [11], as well as the validation method based on the clustered cloud model [12], are only applicable to static outputs. The Bayesian inference-based method focuses on assessing the credibility at the system level [7]. The method based on transfer learning and discrete sequence generative adversarial networks could be used to validate the small-sample dynamic outputs [13]. However, it is hard to meet the data dependency required by transfer learning in practice. There is still a lack of effective validation methods for small-sample dynamic outputs.
Validation methods for dynamic outputs have evolved from directly measuring the original data to approaches based on data features or even data models. Common methods of measuring original data include grey relational analysis (GRA), and theil inequality coefficient (TIC) [14]. However, directly measuring the consistency of original data often reveals only superficial similarities and fails to capture the intrinsic characteristics of the data. Additionally, this approach can be significantly influenced by the uncertainty in the data. Validation methods based on data features assess the consistency of the features such as mean, standard deviation, kurtosis, period, autocorrelation coefficient, and spectral density [9], which are designed and selected according to the specific tasks and types of data involved. Nevertheless, it is challenging to propose guidelines for selecting features for specific systems. Model-based validation methods, in contrast, compare the intrinsic structure and patterns of the data, providing a deeper level of consistency information [15]. These methods can effectively learn the intrinsic structure and capture underlying trends, even in small samples. Therefore, this paper aims to explore a model-based validation method for small-sample dynamic outputs.
The representative models used for time-series modeling and validation mainly include the autoregressive (AR) model [16], the vector AR moving average (VARMA) model [17], the Markov chain model [18], the hidden Markov model [19], and the Gaussian process regression (GPR) model [20]. GPR exhibits strong performance in modeling data uncertainty and complex nonlinear relationships [21,22]. By employing various kernel functions, it can effectively analyze different types of time series data [23]. Thus, we model both the simulation outputs and the reference data by GPR and analyze the consistency between the GPR models.
This paper proposes a validation method based on GPR for small-sample dynamic outputs. The remainder of this paper is structured as follows: Section 2 primarily describes and analyzes the problem of validating small-sample dynamic outputs. Section 3 elaborates a small-sample validation method based on GPR. Within this section, Subsection 3.1 establishes a validation framework based on Bayes statistics, outlining a two-stage validation process that includes segmented GPR and consistency measurement. These two stages are thoroughly described in Subsection 3.2 and Subsection 3.3, respectively. Section 4 provides two examples to demonstrate the effectiveness of the proposed method in describing uncertainty and validating small samples. Finally, the conclusions and future work are presented in Section 5.
The outputs of real systems and simulation models often exhibit uncertainty due to variations in input parameters and the intrinsic randomness of the system. In such cases, it is inaccurate to validate uncertain simulation results based on the outputs of a single experiment. To address this issue, multiple experiments are typically conducted, allowing for the validation of simulation results that include multi-sample simulation outputs and multi-sample reference data. However, complex simulation models cannot be repeated frequently due to high computational costs and extensive time consumption. Moreover, the number of experiments that can be conducted on many real systems is even more limited, considering the costs in terms of money and time. This leads to a common issue: the validation data are often small samples. Hence, our aim is to validate univariate, dynamic, and small-sample outputs.
Suppose that R denotes the real system while S denotes the simulation model. Let R and S run $ {N}_{R} $ times and $ {N}_{S} $ $ \left({N}_{R},{N}_{S}\leq 30\right) $ times respectively. $ {\boldsymbol{Y}}_{R}=\left[\boldsymbol{y}_{R}^{1},\boldsymbol{y}_{R}^{2},\cdots ,\boldsymbol{y}_{R}^{{N}_{R}}\right] $ is an $ {N}_{R}\times M $ matrix which represents $ {N}_{R} $ samples of reference data across $ M $ timestamps. Similarly, $ {\boldsymbol{Y}}_{S}= \left[\boldsymbol{y}_{S}^{1},\boldsymbol{y}_{S}^{2},\cdots ,\boldsymbol{y}_{S}^{{N}_{S}}\right] $ represents $ {N}_{S} $ samples of simulation outputs across $ M $ timestamps. A single sample $ {\boldsymbol{y}}^{n} $ $ \left(n=1,2,\cdots ,{N}_{R}+{N}_{S}\right) $ within $ {\boldsymbol{Y}}_{R} $ or $ {\boldsymbol{Y}}_{S} $ can be denoted as $ {\boldsymbol{y}}^{n}=\left[y_{{t}_{1}}^{n},y_{{t}_{2}}^{n},\cdots ,y_{{t}_{M}}^{n}\right] $ where $ {t}_{1},{t}_{2},\cdots ,{t}_{M} $ are the time points.
Simulation result validation methods obtain the simulation credibility by measuring the consistency between the simulation outputs and the reference data. Suppose that $ C\left({\boldsymbol{Y}}_{S},{\boldsymbol{Y}}_{R}\right) $ denotes the simulation credibility, and let $ C\left({\boldsymbol{Y}}_{S},{\boldsymbol{Y}}_{R}\right)\in [0,1] $. $ C\left({\boldsymbol{Y}}_{S},{\boldsymbol{Y}}_{R}\right)=1 $ when $ {\boldsymbol{Y}}_{S} $ and $ {\boldsymbol{Y}}_{R} $ are completely consistent, and $ C\left({\boldsymbol{Y}}_{S},{\boldsymbol{Y}}_{R}\right)\rightarrow 0 $ as the consistency between $ {\boldsymbol{Y}}_{S} $ and $ {\boldsymbol{Y}}_{R} $ diminishes. Since model-based validation methods capture the intrinsic characteristics of the validation data and require only small-sample data [24], we intend to explore a model-based validation method for small-sample dynamic outputs.
Subsequently, we will study a model-based validation method for small-sample dynamic outputs from two aspects:
(i) We will select an appropriate regression method to model the small-sample dynamic outputs and to capture both the intrinsic characteristics and the uncertainty in the validation data.
GPR is not only applicable to modeling complex data, such as nonlinear and non-smooth data, but also capable of describing the uncertainty in data [25,26]. Hence, the relationship between the outputs and their influences can be fitted by GPR to capture the intrinsic structure and patterns of the validation data. Assuming that factors other than time do not influence the outputs, and given two sets of validation data and the corresponding timestamps, we employ two distinct GPR models to individually describe the reference data $ {\boldsymbol{Y}}_{R} $ and the simulation outputs $ {\boldsymbol{Y}}_{S} $. The elements in $ {\boldsymbol{Y}}_{R} $ and $ {\boldsymbol{Y}}_{S} $ can be respectively regarded as $ {N}_{R}\times M $ samples and $ {N}_{S}\times M $ samples obtained from their corresponding models.
(ii) We will measure the consistency between the simulation outputs and the reference data based on the characteristics captured by the GPR models.
The posterior distributions obtained by GPR describe both the intrinsic characteristics and the uncertainty in the data. According to Bayesian statistical method [27,28], we use the posterior distributions to make inferences about the consistency between two sets of validation data.
In response to the above two aspects, we propose a small-sample validation framework based on Bayesian statistics, as shown in Fig. 1, which mainly consists of two parts: segmented GPR and consistency measurement.
We adopt GPR models, suitable for modeling small-sample time series, to model the validation data. Given the same prior distribution, we model the simulation outputs and the reference data separately by segmented GPR and obtain the posterior distributions.
The consistency between the simulation outputs and the reference data is assessed by comparing the consistency of their posterior distributions in two ways. Firstly, we compare the central tendencies of the two sets of data support by measuring the consistency of the posterior mean vectors. Secondly, we compare the distribution ranges of the two sets of data exhibit by measuring the consistency of the prediction intervals.
A Gaussian process (GP) is a stochastic process composed of Gaussian random variables, where any finite number of these variables follow a joint Gaussian distribution, as shown in Fig. 2.
According to the central limit theorem, the superposition of a large number of independent random variables exhibits a normal distribution when influenced by multiple small, independent uncertainty factors [29]. Therefore, the output data obtained from a single run of a simulation model or a real system can usually be regarded as a single realization of a GP. When the response of a real system or a simulation model is a single dynamic variable with uncertainty, the outputs obtained by repeatedly running are multi-sample time series with the data at each time point following a Gaussian distribution. GPR is a regression algorithm based on a Bayesian approach, relying on observations to model and predict data through the definition of a GP prior [30]. Hence, we could train two independent GPR models to model the simulation outputs and the reference data.
A GP is specified by a mean function $ \mu \left(t\right) $ and a covariance function $ k\left(t,t'\right) $, thereby the validation data are modeled as follows:
$ \begin{cases} y_{t}^{n}=f\left(t\right)+\varepsilon \\f\sim \mathrm{GP}\left(\mu \left(t\right),k\left(t,t'\right)\right)\end{cases} $
where $ y_{t}^{n} $ is the element in $ {\boldsymbol{y}}^{n} $ corresponding to the time point $ t $; $ t $ and $ t' $ are any two points in $ \left[{t}_{1},{t}_{2},\cdots ,{t}_{M}\right] $; $ f $ is the unknown function which models the potential relationship between the validation data and time; $ \varepsilon \sim {\mathrm{N}}\left(0,{\sigma }^{2}\right) $ is the variance, assumed to be independently and identically distributed at different points, indicating the uncertainty in the validation data. The mean function $ \mu \left(t\right) $ defines the expected value of the function given the input $ t $ and is often assumed to be constant or a linear function. The covariance function $ k\left(t,t'\right) $ (also known as the kernel function) defines the covariance between the function values at any two time points $ t $ and $ t' $, capturing the dependencies and structure in the data.
GPR infers the posterior distribution from the prior distribution and the observations. The mean vector and covariance matrix of the posterior distribution can be considered as the characteristics of the time series. The mean vector provides the optimal estimate for each time point and represents the central tendency of the time series. The covariance matrix describes the correlation between the time points and reflects the uncertainty in the time series.
To model the small-sample validation data based on the one-dimensional GPR model, the data need to be preprocessed firstly. Subsequently, we apply the same prior distribution to the unknown functions, selecting an appropriate kernel function based on the characteristics of the validation data. Next, we train and optimize two distinct GPR models, adjusting the model parameters accordingly. After that, two sets of posterior distributions are obtained, which are then used for further consistency analysis. The process is shown in Fig. 3.
Data preprocessing includes three parts: matrix transformation, normalization and segmentation. Since we model the validation data based on the one-dimensional GPR, before being input into the regression model, $ {\boldsymbol{Y}}_{R} $ needs to be converted from an $ {N}_{R}\times M $ matrix to a vector of length $ {N}_{R}M $, denoted as $ {\boldsymbol{Y}}_{R}' $. Similarly, $ {\boldsymbol{Y}}_{S} $ needs to be converted from an $ {N}_{S}\times M $ matrix to a vector of length $ {N}_{S}M $, denoted as $ {\boldsymbol{Y}}_{S}' $. Two sets of validation data and their corresponding time data, $ {\boldsymbol{T}}_{R}' $ and $ {\boldsymbol{T}}_{S}' $ are as follows:
$ \left\{\begin{aligned}&{\boldsymbol{Y}}_{R}'={\left[y_{R{t}_{1}}^{1},y_{R{t}_{1}}^{2},\cdots ,y_{R{t}_{1}}^{{N}_{R}},\cdots ,y_{R{t}_{M}}^{1},y_{R{t}_{M}}^{2},\cdots ,y_{R{t}_{M}}^{{N}_{R}}\right]}_{\left(1\times {N}_{R}M\right)}\\&{\boldsymbol{T}}_{R}'={\left[{t}_{1},\cdots ,{t}_{1},\cdots ,{t}_{M},\cdots ,{t}_{M}\right]}_{\left(1\times {N}_{R}M\right)}\end{aligned} \right., $
$ \left\{\begin{aligned}&{\boldsymbol{Y}}_{S}'={\left[y_{S{t}_{1}}^{1},\cdots ,y_{S{t}_{1}}^{{N}_{S}},\cdots ,y_{S{t}_{M}}^{1},\cdots ,y_{S{t}_{M}}^{{N}_{S}}\right]}_{\left(1\times {N}_{S}M\right)}\\&{\boldsymbol{T}}_{S}'={\left[{t}_{1},\cdots ,{t}_{1},\cdots ,{t}_{M},\cdots ,{t}_{M}\right]}_{\left(1\times {N}_{S}M\right)}\end{aligned}\right., $
where $ y_{Rt}^{{n}_{1}}\left({n}_{1}=1,2,\cdots ,{N}_{R},t={t}_{1},{t}_{2},\cdots ,{t}_{M}\right) $ is the element in the $ {n}_{1}\text{th} $ sample of reference data corresponding to the time point $ t $; and $ y_{St}^{{n}_{2}}\left({n}_{2}=1,2,\cdots ,{N}_{S}\right) $ is the element in the $ {n}_{\text{2}}{\mathrm{th}} $ sample of simulation outputs corresponding to the time point $ t $.
Additionally, in order to improve model performance and stability, normalization is performed on $ {\boldsymbol{Y}}_{R}' $ and $ {\boldsymbol{Y}}_{S}' $ as follows:
$ y_{R{\mathrm{\text{norm}}}}^{i}=\frac{y_{R}^{i}-\min \left\{y_{R}^{i},y_{S}^{j}\right\}}{\max \left\{y_{R}^{i},y_{S}^{j}\right\}-\min \left\{y_{R}^{i},y_{S}^{j}\right\}} , $
$ y_{S{\mathrm{\text{norm}}}}^{j}=\frac{y_{S}^{j}-\min \left\{y_{R}^{i},y_{S}^{j}\right\}}{\max \left\{y_{R}^{\mathrm{i}},y_{S}^{j}\right\}-\min \left\{y_{R}^{i},y_{S}^{j}\right\}}, $
where $ y_{R{\mathrm{\text{norm}}}}^{i} $ and $ y_{S{\mathrm{\text{norm}}}}^{i} $ are the elements in the normalized validation data $ {\boldsymbol{Y}}_{R{\mathrm{\text{norm}}}} $ and $ {\boldsymbol{Y}}_{S{\mathrm{\text{norm}}}} $; $ i=1,2,\cdots ,{N}_{R}M $; $ j=1,2,\cdots ,{N}_{S}M $; $ y_{R}^{i} $ and $ y_{S}^{j} $ are the elements in the validation data $ {\boldsymbol{Y}}_{R}' $ and $ {\boldsymbol{Y}}_{S}' $; $ \min \left\{y_{R}^{i},y_{S}^{j}\right\} $ denotes the minimum element in $ {\boldsymbol{Y}}_{R}' $ and $ {\boldsymbol{Y}}_{S}' $; and $ \max \left\{y_{R}^{i},y_{S}^{j}\right\} $ denotes the maximum element in $ {\boldsymbol{Y}}_{R}' $ and $ {\boldsymbol{Y}}_{S}' $.
Due to the numerous matrix inverse operations involved in the training process of the GPM, the calculations would increase with the size of the dataset. The complexity of GPR is $ O\left({l}^{3}\right) $ when the length of the sequence is $ l $ [31]. The validation data in simulation validation problem are usually long time series which would result in a large amount of computation and long time-consumption. To solve this problem, we perform segmented GPR with a sliding window for data segmentation.
The size of the sliding window mainly depends on the characteristics of the data itself. If the time series is relatively smooth or exhibits no obvious trend changes, the window can be set to a fixed size based on the amount of data. For example, for periodic data, the size of the window can be determined according to the period of the data. If the time series shows an obvious trend change, the data can be divided according to the trend changes at different stages using windows of varying sizes. Positions in the time series where there is a clear change in trend or a jump in value are reasonable segmentation points, and the window size can be adjusted accordingly. Additionally, the number of samples in each segment should not be too small; otherwise, the fitting effect of GPR may be suboptimal. On the other hand, the number of samples in each segment should not be too large, as it may affect computational efficiency. Typically, the number of samples per segment can range between 50 and 200. This ensures both a good fitting effect and manageable computational cost.
Assuming that the length of the window is $ w $ timestamps, we slide $ w $ timestamps each time so that two sets of normalized validation data, $ {\boldsymbol{Y}}_{R{\mathrm{\text{norm}}}} $ and $ {\boldsymbol{Y}}_{S{\mathrm{\text{norm}}}} $, are both divided into $ q $ segments. The $ {x}{{\mathrm{th}}} $ segment of the validation data, $ {\boldsymbol{Y}}_{Rx} $ and $ {\boldsymbol{Y}}_{Sx} $, and their corresponding time data, $ {\boldsymbol{T}}_{Rx} $ and $ {\boldsymbol{T}}_{Sx} $, are as follows:
$ \left\{\begin{aligned}&{\boldsymbol{Y}}_{Rx}={\left[y_{R{\mathrm{\text{norm}}}}^{xw+1},\cdots ,y_{R{\mathrm{\text{norm}}}}^{xw+{N}_{R}w}\right]}_{\left(1\times {N}_{R}w\right)}\\&{\boldsymbol{T}}_{Rx}={\left[{t}_{xw+1},\cdots ,{t}_{xw+1},\cdots ,{t}_{xw+w},\cdots ,{t}_{xw+w}\right]}_{\left(1\times {N}_{R}w\right)}\end{aligned}\right. , $
$ \left\{\begin{aligned}&{\boldsymbol{Y}}_{Sx}={\left[y_{S{\mathrm{\text{norm}}}}^{xw+1},\cdots ,y_{S{\mathrm{\text{norm}}}}^{xw+{N}_{S}w}\right]}_{\left(1\times {N}_{S}w\right)}\\&{\boldsymbol{T}}_{Sx}={\left[{t}_{xw+1},\cdots ,{t}_{xw+1},\cdots ,{t}_{xw+w},\cdots ,{t}_{xw+w}\right]}_{\left(1\times {N}_{S}w\right)}\end{aligned}\right., $
where $ y_{R{\mathrm{\text{norm}}}}^{i}\left(xw+1\leq i\leq xw+{N}_{R}w,i\in \mathbb{N}\right) $ is the $ {i}{\text{th}} $ element in $ {\boldsymbol{Y}}_{R\text{norm}} $; $ x=1,2,\cdots ,q $, $ q=\dfrac{M}{w} $; and $ y_{S\text{norm}}^{j}\left(xw+1\leq j\leq xw+{N}_{S}w,\right. $ $ \left.j\in {{\bf{N}}}\right) $ is the $ {j}{\text{th}} $ element in $ {\boldsymbol{Y}}_{S\text{norm}} $.
Let the unknown functions corresponding to $ \left({\boldsymbol{T}}_{Rx},{\boldsymbol{Y}}_{Rx}\right) $ and $ \left({\boldsymbol{T}}_{Sx},{\boldsymbol{Y}}_{Sx}\right) $ be denoted as $ {f}_{Rx} $ and $ {f}_{Sx} $, then we aim to learn the posterior distributions of them. Due to the same methods and processes for modeling each segment of the validation data, including the reference data and the simulation outputs, the following is a unified description. The $ {x}{{\mathrm{th}}} $ segment of the validation data, the time data and the function are respectively denoted as $ {\boldsymbol{Y}}_{x} $, $ {\boldsymbol{T}}_{x} $ and $ {f}_{x} $, and $ N $ is the sample size of the validation data.
We need to define the prior distribution of the function after data preprocessing. Let the prior distribution of the unknown function $ {f}_{x} $ be a GP with zero mean function and a kernel function, denoted as $ {f}_{x}\sim \mathrm{GP}\left(0,{k}_{x}\left(t,t'\right)\right) $. GPR effectively learns complex trends and nonlinear relationships in data by utilizing appropriate kernel functions, making it adaptable to a wide variety of data forms. Through the selection of suitable kernel functions, GPR can capture intricate relationships in time series data, including periodicity, trends, and self-similarity. Common kernel functions include the Gaussian kernel (also known as the radial basis function, RBF), Matern kernel, linear kernel, polynomial kernel, and periodic kernel. The Gaussian kernel is particularly well-suited for modeling locally smooth and nonlinear time series data:
$ {k}_{\text{RBF}}\left(t,t'\right)={\sigma }_{1}^{2}\exp \left(-\frac{{\left|\left|t-t'\right|\right|}^{2}}{2{l}_{1}^{2}}\right) $
where $ {\sigma }_{1}^{2} $ is the variance, controlling the vertical variation of the function; and $ {l}_{1} $ is the length scale parameter, controlling the horizontal smoothness of the function.
And it often satisfies the fitting requirements of validation data. However, when the validation data exhibit significant fluctuations or abrupt changes, the Matern kernel is more appropriate:
${k}_{\text{Matern}}\left(t,t'\right)= \frac{{2}^{1-\nu }}{\mathit{\Gamma }(\nu )}{\left(\frac{\sqrt{2\nu }\left| t-t'\right| }{{l}_{2}}\right)}^{\nu }{K}_{\nu }\left(\frac{\sqrt{2\nu }\left| t-t'\right| }{{l}_{2}}\right) $
where $ \nu $ is the smoothness parameter; $ \mathit{\Gamma }(\nu ) $ is the Gamma function; $ {l}_{2} $ is the length scale; and $ {K}_{\nu } $ is the modified Bessel function.
In addition, other kernel functions, such as the linear kernel, polynomial kernel, and periodic kernel, can be chosen based on prior knowledge or observed trends in the data.
The prior distribution is used to initialize the GPR models. Then we can optimize the GPR model with the preprocessed validation data $ {\boldsymbol{Y}}_{x} $ and the corresponding time data $ {\boldsymbol{T}}_{x} $ as the training data. To ensure the GPR models accurately fits the validation data, the model’s parameters, which include the hyperparameters in the kernel function $ {k}_{x}\left(t,t'\right) $ and the variance $ {\sigma }^{2} $, are meticulously optimized by the maximum likelihood method. This optimization process involves adjusting the parameters through a gradient descent algorithm, iteratively refining them until the parameters that maximize the log marginal likelihood (LML) are identified. This method improves the reliability and effectiveness of the GPR models in capturing the essential characteristics of the validation data. LML is calculated as follows:
$ \begin{gathered}[b]\lg p\left({\boldsymbol{Y}}_{x}|{\boldsymbol{T}}_{x},{\theta }_{x}\right)=-\frac{1}{2}\left({\boldsymbol{Y}}_{x}^{\mathrm{T}}{\left({\boldsymbol{K}}_{x}+{\sigma }^{2}\boldsymbol{I}\right)}^{-1}{\boldsymbol{Y}}_{x}\right.+\\ \left.\lg \left| {\boldsymbol{K}}_{x}+{\sigma }^{2}\boldsymbol{I}\right| +wN\lg \left(2\text{π} \right)\right)\end{gathered} $
where $ {\boldsymbol{Y}}_{x} $ and $ {\boldsymbol{T}}_{x} $ are the output and input of the GPR model; $ {\theta }_{x} $ denotes the parameters to be optimized; $ {\boldsymbol{K}}_{x} $ is the covariance matrix of the training time points $ \boldsymbol{T}_{x}^{} $, consisting of the elements calculated by $ {k}_{x}\left(t,t'\right),t,t'\in {\boldsymbol{T}}_{x} $; and $ {\sigma }^{2} $ is the variance.
After that, we can input test time points into the optimized model to obtain the posterior distribution. Due to our desire to capture the essential characteristics of the validation data, we make the test time points $ \boldsymbol{T}_{x}^{*} $ consist of the non-repetitive timestamps corresponding to the validation data, as follows:
$ \boldsymbol{T}_{x}^{*}={\left[t_{xw+1}^{},t_{xw+2}^{},\cdots ,t_{xw+w}^{}\right]}_{\left(1\times w\right)} .$
The posterior distribution of the outputs, denoted as $ \boldsymbol{y}_{x}^{*} $, follows a Gaussian distribution with posterior mean $ \boldsymbol{\mu }_{x}^{*} $ and posterior covariance $ {\boldsymbol{\varSigma }}_{x}^{*} $:
$ \boldsymbol{y}_{x}^{*}|\boldsymbol{T}_{x}^{*},{\boldsymbol{T}}_{x},{\boldsymbol{Y}}_{x}\sim {\mathrm{N}}\left(\boldsymbol{\mu }_{x}^{*},{\boldsymbol{\varSigma }}_{x}^{*}\right) , $
$ \boldsymbol{\mu }_{x}^{*}=\boldsymbol{K}_{x}^{*}{\left({\boldsymbol{K}}_{x}+{\sigma }^{2}\boldsymbol{I}\right)}^{-1}{\boldsymbol{Y}}_{x} , $
$ {\boldsymbol{\varSigma}} _{x}^{*}=\boldsymbol{K}_{x}^{**}-\boldsymbol{K}_{x}^{*}{\left({\boldsymbol{K}}_{x}+{\sigma }^{2}\boldsymbol{I}\right)}^{-1}{\left(\boldsymbol{K}_{x}^{*}\right)}^{\text{T}}+{\sigma }^{2}\boldsymbol{I}, $
where $ \boldsymbol{K}_{x}^{*} $ is the covariance between the test time points $ \boldsymbol{T}_{x}^{*} $ and the training time points $ \boldsymbol{T}_{x}^{} $; and $ \boldsymbol{K}_{x}^{**} $ is the covariance among the test time points $ \boldsymbol{T}_{x}^{*} $.
Following Bayes statistics method, we place the same prior distribution $ {\mathrm{GP}}(0,{k}_{x}\left(t,t'\right)) $ on the unknown functions $ f_{Sx} $ and $ f_{Rx} $; then, given the $ {x}{{\mathrm{th}}} $ segment of the validation data and their time data $ \left({\boldsymbol{T}}_{Sx},{\boldsymbol{Y}}_{Sx}\right) $ and $ \left({\boldsymbol{T}}_{Rx},{\boldsymbol{Y}}_{Rx}\right) $, we obtain their posterior distributions:
$ \boldsymbol{y}_{Rx}^{*}\sim {\mathrm{N}}(\boldsymbol{\mu }_{Rx}^{*},{\boldsymbol{\varSigma }}_{Rx}^{*}) , $
$ \boldsymbol{y}_{Sx}^{*}\sim {\mathrm{N}}(\boldsymbol{\mu }_{Sx}^{*},{\boldsymbol{\varSigma }}_{Sx}^{*}) . $
Next, we make inferences about the consistency between $ {\boldsymbol{Y}}_{Sx} $ and $ {\boldsymbol{Y}}_{Rx} $ by comparing the posterior distributions; finally, the validation results for each segment of data are synthesized to indicate the simulation model’s credibility.
The posterior mean vector represents the central tendency of the validation data, while the prediction interval reflects the range of data distribution, taking into account the uncertainty of the validation data. If two posterior distributions are similar in terms of central tendency and distribution range, it can be inferred that both the intrinsic characteristics and the uncertainty in two sets of validation data are consistent. Thus, we simplify the task of comparing the consistency of posterior distributions to quantifying the consistency of mean vectors and prediction intervals.
We compare the central tendencies of the function behavior corresponding to the simulation outputs and the reference data by comparing the posterior mean vectors. The Euclidean distance and cosine distance can respectively measure differences in position and direction between two vectors, respectively. Hence, to assess the consistency of the two posterior mean vectors $ \boldsymbol{\mu }_{Sx}^{*} $ and $ \boldsymbol{\mu }_{Rx}^{*} $, we comprehensively measure the differences between them by Euclidean distance and cosine distance, as follows:
$ {D}_{\mathrm{eu}-x}={\left|\left|\boldsymbol{\mu }_{Sx}^{*}-\boldsymbol{\mu }_{Rx}^{*}\right|\right|}_{2} , $
$ {D}_{\cos -x}=\frac{\boldsymbol{\mu }_{Sx}^{*}\cdot \boldsymbol{\mu }_{Rx}^{*}}{\left|\left|\boldsymbol{\mu }_{Sx}^{*}\right|\right|\left|\left|\boldsymbol{\mu }_{Rx}^{*}\right|\right|} . $
We aim to map the range of consistency between the two posterior mean vectors to [0,1]. The smaller the $ {D}_{\mathrm{eu}-x} $, the higher the degree of consistency. The larger the $ {D}_{\cos -x} $, the higher the degree of consistency. Considering that $ {D}_{\mathrm{eu}-x}\in [0,+\mathrm{\infty }] $ and $ {D}_{\cos -x}\in [-1,1] $, we convert the distance metrics as follows:
$ {C}_{mx}=\frac{1}{2}\left(\frac{1}{1+{D}_{\mathrm{eu}-x}}+\frac{1+{D}_{\cos -x}}{2}\right) $
where $ {C}_{mx} $ denotes the consistency of the posterior mean vectors.
We compare the coverage range of two prediction intervals to quantify and compare the uncertainty in the simulation outputs and the reference data. For computational simplicity and intuitive understanding, we consider the area of the prediction interval as the sum of the widths of the intervals at each point. Furthermore, we use the ratio of the intersection area and union area to indicate the consistency of the prediction intervals $ {C}_{{\mathrm{PI}}x} $, as follows:
$ {C}_{{\mathrm{PI}}x}=\frac{{{\mathrm{Area}}}_{Ix}}{{{\mathrm{Area}}}_{Ux}} , $
$ {{\mathrm{Area}}}_{Ix}=\sum \max (\min ({U}_{Sd},{U}_{Rd})-\max ({L}_{Sd},{L}_{Rd}),0) , $
$ {{\mathrm{Area}}}_{Ux}=\sum \left(\left({U}_{Sd}-{L}_{Sd}\right)+\left({U}_{Rd}-{L}_{Rd}\right)\right)-{{\mathrm{Area}}}_{Ix}, $
where $ {{\mathrm{Area}}}_{Ix} $ denotes the intersection area; $ {{\mathrm{Area}}}_{Ux} $ denotes the union area; $ [{L}_{Sd},{U}_{Sd}] $ and $ [{L}_{Rd},{U}_{Rd}] $ are the prediction intervals at time point $ t_{d}^{}(d=xw+1,xw+2,\cdots , xw+w) $ corresponding to the simulation outputs and the reference data:
$ \left\{\begin{aligned}{L}_{Sd}&={\mu }_{Sd}-{z}_{\alpha /2}\cdot {\sigma }_{Sd}\\{U}_{Sd}&={\mu }_{Sd}+{z}_{\alpha /2}\cdot {\sigma }_{Sd}\end{aligned}\right. , $
$ \left\{\begin{aligned}{L}_{Rd}&={\mu }_{Rd}-{z}_{\alpha /2}\cdot {\sigma }_{Rd}\\{U}_{Rd}&={\mu }_{Rd}+{z}_{\alpha /2}\cdot {\sigma }_{Rd}\end{aligned}\right., $
where $ {z}_{\alpha /2} $ is the z-score corresponding to the $ \alpha /2 $ quantile of the standard normal distribution, associated with a confidence level $ \left(1-\alpha \right)\times 100\% $; $ {\mu }_{Sd} $ and $ {\mu }_{Rd} $ are the posterior means at $ t_{d}^{} $, which are derived from the posterior mean vectors $ \boldsymbol{\mu }_{Sx}^{*} $ and $ \boldsymbol{\mu }_{Rx}^{*} $; $ {\sigma }_{Sd} $ and $ {\sigma }_{Rd} $ are the posterior standard deviations at $ t_{d}^{} $, which are computed as the square roots of the diagonal elements of the posterior covariance matrices $ {\boldsymbol{\varSigma }}_{Sx}^{*} $ and $ {\boldsymbol{\varSigma }}_{Rx}^{*} $, respectively.
Then, the validation result of the $ {x}{{\mathrm{th}}} $ segment of the validation data, denoted as $ C({\boldsymbol{Y}}_{Sx},{\boldsymbol{Y}}_{Rx}) $, can be calculated by
$ C({\boldsymbol{Y}}_{Sx},{\boldsymbol{Y}}_{Rx})={\omega }_{1}{C}_{mx}+{\omega }_{2}{C}_{{\mathrm{PI}}x} $
where $ {\omega }_{1} $ and $ {\omega }_{2} $ control the contribution of two parts of consistency to the validation result; $ {\omega }_{1}+{\omega }_{2}=1,0 \lt {\omega }_{1}, {\omega }_{2} \lt 1 $. In general, assessors prioritize the central tendency of the validation data over its uncertainty. Therefore, $ {\omega }_{1} $ is typically assigned a larger value than $ {\omega }_{2} $. The higher the consistency requirement for uncertainty, the larger $ {\omega }_{2} $.
Since we have completed segmented model validation, we synthesize the validation result of each segment $ C({\boldsymbol{Y}}_{Sx},{\boldsymbol{Y}}_{Rx}) $ $ (x=1,2,\cdots ,q) $ into the final validation result $ C({\boldsymbol{Y}}_{S},{\boldsymbol{Y}}_{R}) $ based on the weighted sum method, as follows:
$ C({\boldsymbol{Y}}_{S},{\boldsymbol{Y}}_{R})=\sum \limits_{x=1}^{q}{\lambda }_{x}C({\boldsymbol{Y}}_{Sx},{\boldsymbol{Y}}_{Rx}) $
where $ {\lambda }_{x}(x=1,2,\cdots ,q) $ denote the weight of $ C({\boldsymbol{Y}}_{Sx},{\boldsymbol{Y}}_{Rx}) $, generally set equal to $ 1/q $. Besides, the weights could be adjusted if we focus more on a certain segment of the outputs.
In this section, we evaluate the effectiveness of the proposed validation method for uncertainty description by a numerical example.
Given two simulation models:
$ \left\{\begin{aligned}{{y}}_{1}&=1+{a}_{1}\sin (0.2\text{π} {x})+{b}_{1}\\{{y}}_{2}&=1+{a}_{2}\sin (0.2\text{π} {x})+{b}_{2}\end{aligned}\right. $
where $ 0\leq x\leq 19,x\in {{\bf{R}}} $. We set the parameters of the two simulation models $ \boldsymbol{\zeta }=[{a}_{1},{a}_{2},{b}_{1},{b}_{2}] $ as
$ \begin{gathered}[b] {\boldsymbol{\zeta }}_{1}\colon {a}_{1},{a}_{2}\sim {\mathrm{N}}(1,0.25), {b}_{1},{b}_{2}\sim {\mathrm{N}}(1,0.25),\\ {\boldsymbol{\zeta }}_{2}\colon {a}_{1}\sim {\mathrm{N}}(1,0.25),{a}_{2}\sim {\mathrm{N}}(1,0.5),{b}_{1},{b}_{2}\sim {\mathrm{N}}(1,0.25),\\ {\boldsymbol{\zeta }}_{3}\colon {a}_{1},{a}_{2}\sim {\mathrm{N}}(1,0.25), {b}_{1}\sim {\mathrm{N}}(1,0.5),{b}_{2}\sim {\mathrm{N}}(1,0.25),\end{gathered} $
and carry out three sets of experiments. The two models are repeatedly run 10 times to get two sets of observations in each experiment. After training their corresponding GPR models, the posterior distributions are obtained, as shown in Fig. 4.
Each parameter in $ {\zeta } $ follows a Gaussian distribution, and the larger the variance the greater the uncertainty in the parameter. As seen in Fig. 4, the prediction interval is larger when the uncertainty in the model parameters is greater. Hence, the GPR model is able to effectively describe the uncertainty in the observations.
In this section, we evaluate the effectiveness of the proposed validation method for small-sample validation by an application example on a force load identification model of a flight vehicle. Considering the strain error, we validate the dynamic outputs under uncertainty, focusing on the axial force $ F $ as the target variable.
Let the model be run 20 times with a random strain error of 1% to obtain 20 validation data samples, with 10 of the samples designated as a set of reference data and the other 10 as a set of simulation outputs. Additionally, let the model be run 10 times with a random strain error of 3% to obtain 10 samples as another set of simulation outputs. The reference data are labeled as “Data_R” and two sets of simulation outputs corresponding to 1% and 3% strain errors are labeled as “Data_S1”and “Data_S2”, respectively. The information about the validation data is shown in Table 1.
We use axial force data from 1.01 s to 2.00 s as the validation data, which are dynamic outputs at intervals of 0.01 s. Each run produces a time series spanning 100 timestamps. Five samples of validation data each with 1% and 3% strain error are randomly selected, and the distribution of the data is shown in Fig. 5.
As shown in Fig. 5, the axial force data exhibits periodic behavior with a period of 0.50 s. Additionally, the larger the strain error, the greater the uncertainty in the model output. Then, we validate the two sets of simulation outputs separately and show the validation process. For simplicity and efficiency, we set the sliding window size at 50 and adopt a Gaussian kernel function to fit the validation data. After training the GPR models, we obtain the posterior distributions. The posterior distributions of the first segment of reference data and two sets of simulation outputs are shown in Fig. 6.
Fig. 6 shows that the range of prediction intervals of “Data_S2” is significantly larger. This again indicates that the prediction intervals effectively reflect the uncertainty in the observations.
After that, we can quantify the consistency of the posterior mean vectors and prediction intervals by (17)−(24). Next, the validation results of each segment of data are calculated by (25) with $ {\omega }_{1}=0.7 $, $ {\omega }_{2}=0.3 $. Furthermore, let $ {\lambda }_{1}={\lambda }_{2}=0.5 $ in (26) and obtain the final validation results, as seen in Table 2.
According to the information in Table 1 and the results shown in Table 2, it is evident that the credibility of the simulation decreases as the divergence in parameter distributions increases. This observation is consistent with theoretical expectations, thereby confirming the effectiveness of the proposed validation method in achieving reasonable results.
To demonstrate that segmented GPR has minimal or no impact on the overall validation results, validations are conducted using “Data_S2” as the simulation outputs in two cases with window sizes of 20 and 50. The corresponding validation results and times are presented in Table 3.
It shows that the overall validation results are very close when setting different window sizes for segmented GPR. Besides, the validation process is less time-consuming with a smaller window. This suggests that segmented GPR has almost no impact on the final model validation results but can significantly enhance the validation efficiency.
Furthermore, we test the performance of the proposed method in validating small samples. Let the force load identification model run 100 times under a strain error of 1%, and the model outputs form the reference data set. Then, let the model run 100 times under a strain error of 3%, and the model outputs form the simulation data set. Following Table 4, in the corresponding validation data set, randomly select a specified number of validation data samples for simulation results validation. The final validation results are also shown in Table 4.
It can be seen that the validation results are very similar, whether only the sample size of simulation outputs is increased or the sample sizes of both simulation outputs and reference data are increased. Therefore, the proposed method can be effectively applied to simulation result validation under small samples.
To verify the superiority of the proposed method, we compare the proposed validation method based on GPR with another simulation dynamic output validation method based on data features proposed in [9]. The validation method based on data features is designed to measure dynamic simulation results under uncertainty based on Bayesian hypothesis testing and it is assumed that the simulation outputs are consistent with the reference data if the validation result exceeds 0.5.
Referring to the validation data information in Table 1, the model outputs under a 1% strain error are used as the reference data. Besides, the model outputs under 1% and 3% strain errors are used as two sets of simulation data, labeled as “Data_S1” and “Data_S2” respectively. The two methods are applied in three cases with sample sizes of 10, 30, and 50 for both reference data and simulation data. The validation results are shown in Fig. 7.
As can be seen from Fig. 7, the proposed method always yields reasonable validation results even with small sample sizes. However, the validation method based on data features can obtain effective validation results when the samples are large but its performance is poor when validating small samples. In summary, the proposed method demonstrates significant superiority in small-sample validation.
To address the challenge of model validation under small samples, this paper proposes a validation method based on GPR to validate small-sample dynamic outputs. We establish a validation framework based on Bayes statistics, shifting the focus from analyzing validation data to analyzing posterior distributions. The validation process is divided into two main stages: segmented GPR and consistency measurement. We respectively model the simulation outputs and the reference data by segmented GPR. Furthermore, the consistency of the validation data is assessed by comparing the posterior distributions of the GPR models, which effectively describe the intrinsic characteristics and the uncertainty in validation data.
We test the effectiveness of the proposed method in terms of both uncertainty description and validation for small samples through a numerical example and an application example. The method has been demonstrated to achieve reasonable results for validating small samples.
However, the proposed method is currently limited to dynamic and univariate outputs. In reality, outputs may be multivariate and correlated. Future work will focus on extending the small-sample validation to accommodate multivariate and correlated outputs.
1
STEVENSON D E. Verification and validation of complex systems. Proc. of the Intelligent Engineering Systems through Artificial Neural Networks, 2002: 159−164.
2
WANG Y N, LI J Q, SUN H B, et al. A survey on VV&A of large-scale simulations. International Journal of Crowd Science, 2019, 3(1): 63–86.
3
SARGENT R G. Verification and validation of simulation models. Journal of Simulation, 2013, 7(1): 12–24.
4
DURST P J, ANDERSON D T, BETHEL C L. A historical review of the development of verification and validation theories for simulation models. International Journal of Modeling, Simulation, and Scientific Computing, 2017, 8(2): 1730001.
5
XI Z M, YANG R J. Reliability analysis with model uncertainty coupling with parameter and experiment uncertainties: a case study of 2014 V&V challenge problem. Journal of Verification, Validation and Uncertainty Quantification, 2015, 1(1): 011005.
6
SHANKAR S, SANKARAN M. Integration of model verification, validation, and calibration for uncertainty quantification in engineering systems. Reliability Engineering and System Safety, 2015, 138: 194–209.
7
LIN S L, LI W, MA P, et al. Structural modelling and bayesian inference for model validation and confidence extrapolation. Journal of Statistical Computation and Simulation, 2020, 90(2): 211–233.
8
PAN Y L, HE Y Z. Trustworthiness assessment method for guidance simulation system based on DS/AHP and gray cloud clustering. Electronic Measurement Technology, 2017, 40(7): 43–47.
9
YANG M, QIAN X C, LI W. Simulation dynamic output validation method based on data feature. Systems Engineering and Electronic Technology, 2016, 38(2): 457–463. (in Chinese)
10
JIA J P, HE X Q, JIN Y J. Statistics. Beijing: Renmin University Press, 2015.(in Chinese)
11
SONG T. Research on model validation method for small sample and service-oriented tools. Harbin: Harbin Institute of Technology, 2018. (in Chinese)
12
WANG J M, WU Y J. Credibility assessment of small sample data based on clustered cloud model. Journal of System Simulation, 2019, 31(7): 1263–1271. (in Chinese)
13
NIE K, LUAN R P. A simulation model validation method based on data enhancement. Command Control and Simulation, 2019, 41(3): 92–96. (in Chinese)
14
LI W, ZHOU Y C, LIN S L, et al. A review of simulation model validation methods. Journal of System Simulation, 2019, 31(7): 1249–1256. (in Chinese)
15
ANALLA M. Model validation through the linear regression fit to actual versus predicted values. Agricultural Systems, 1998, 57(1): 115–119.
16
MAHARAJ E A. Cluster of time series. Journal of Classification, 2000, 17(2): 297–314.
17
MAHARAJ E A. Comparison and classification of stationary multivariate time series. Pattern Recognition, 1999, 32(7): 1129–1138.
18
RAMONI M, SEBASTIANI P, COHEN P. Bayesian clustering by dynamics. Machine Learning, 2002, 47(1): 91–121.
19
BICEGO M, MURINO V, FIGUEIREDO M A T. Similarity-based clustering of sequences using hidden Markov models. Proc. of the International Conference on Machine Learning and Data Mining in Pattern Recognition, 2003: 86−95.
20
SEEGER M. Gaussian processes for machine learning. International Journal of Neural Systems, 2004, 14(2): 69–106.
21
ROBERTS S, OSBORNE M, EBDEN M, et al. Gaussian processes for time-series modelling. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 2013, 371(1984): 20110550.
22
AIGRAIN S, FOREMAN-MACKEY D. Gaussian process regression for astronomical time series. Annual Review of Astronomy and Astrophysics, 2023, 61(1): 329–371.
23
PALAR P S, PARUSSINI L, BREGANT L, et al. On kernel functions for bi-fidelity Gaussian process regressions. Structural and Multidisciplinary Optimization, 2023, 66(2): 1–22.
24
LI S. Research on the validation method of missile system simulation model. Changsha: National University of Defense Technology, 2003. (in Chinese)
25
MANFREDI P. Conservative Gaussian process models for uncertainty quantification and Bayesian optimization in signal integrity applications. IEEE Trans. on Components Packaging and Manufacturing Technology, 2024, 14(7): 1261–1272.
26
KARVONEN T, WYNNE G, TRONARP F, et al. Maximum likelihood estimation and uncertainty quantification for gaussian process approximation of deterministic functions. SIAM/ASA Journal on Uncertainty Quantification, 2020, 8(3): 926–958.
27
PARK I, AMARCHINTA H K, GRANDHI R V. A Bayesian approach for quantification of model uncertainty. Reliability Engineering and System Safety, 2010, 95(7): 777–785.
28
REBBA R, MAHADEVAN S, HUANG S. Validation and error estimation of computational models. Reliability Engineering and System Safety, 2006, 91(10): 1390–1397.
29
SMITH J. Statistical modeling principles. Journal of Statistical Methods, 2023, 45(2): 234–256.
30
ANDERSON R, BROWN T. Gaussian process regression for dynamic systems. Computational Bayesian Models, 2023, 29(3): 98–121.
31
LIU H, ONG Y S, SHEN X, et al. When Gaussian process meets big data: a review of scalable GPS. IEEE Trans. on Neural Networks and Learning Systems, 2020, 31(11): 4405–4423.
Year 2026 volume 37 Issue 3
PDF
81
48
Cite this Article
BibTeX
Article Info
doi: 10.23919/JSEE.2026.000113
  • Receive Date:2024-05-29
  • Online Date:2026-08-14
  • Published:2026-06-18
Article Data
Affiliations
History
  • Received:2024-05-29
Affiliations
    1Control and Simulation Center, Harbin Institute of Technology, Harbin 150080, China
    2National Key Laboratory of Complex System Modeling and Simulation, Harbin 150080, China

Corresponding:

LI Wei
References
Share
https://castjournals.cast.org.cn/joweb/jsee/EN/10.23919/JSEE.2026.000113
Share to
QR

Scan QR to access full text

Cite this article
BibTeX
Citations
表12种不同金属材料的力学参数

Family
属数
Number of
genus
种数
Number of
species
占总种数比例
Percentage of
total species (%)

Genus
种数
Number of
species
占总种数比例
Percentage of total
species (%)
鹅膏菌科Amanitaceae 2 11 5.26 鹅膏菌属 Amanita 10 4.78
小菇科 Mycenaceae 2 12 5.74 丝盖伞属 Inocybe 5 2.39
多孔菌科 Polyporaceae 8 14 6.70 蜡蘑属 Laccaria 5 2.39
红菇科 Russulaceae 3 23 11.00 小皮伞属 Marasmius 6 2.87
小菇属 Mycena 11 5.26
光柄菇属 Pluteus 5 2.39
红菇属 Russula 17 8.13
栓菌属 Trametes 5 2.39
关闭全屏
  • BibTeX
  • EndNote
  • RefWorks
  • TxT