Abstract

While most calibration methods focus on inferring a set of model parameters that are unknown but assumed to be constant, many models have parameters that have a functional relation with the controllable input variables. Formulating a low-dimensional approximation of these calibration functions allows modelers to use low-fidelity models to explore phenomena at lengths and time scales unattainable with their high-fidelity sources. While functional calibration methods are available for low-dimensional problems (e.g., one to three unknown calibration functions), exploring high-dimensional spaces of unknown calibration functions (e.g., more than ten) is still a challenging task due to its computational cost and the risk for identifiability issues. To address this challenge, we introduce a semiparametric calibration method that uses an approximate Bayesian computation scheme to quantify the uncertainty in the unknown calibration functions and uses this insight to identify what functions can be replaced with low-dimensional approximations. Through a test problem and a coarse-grained model of an epoxy resin, we demonstrate that the introduced method enables the identification of a low-dimensional set of calibration functions with a limited compromise in calibration accuracy. The novelty of the presented method is the ability to synthesize domain knowledge from various sources (i.e., physical experiments, simulation models, and expert insight) to enable high-dimensional functional calibration without the need for prior knowledge on the class of unknown calibration functions.

1 Introduction

Simulation models provide a low-fidelity approximation of physical systems but need to be thoroughly validated before being used to inform design decisions. For example, a high-fidelity response model fh(x) as given by the red curve in Fig. 1 can be approximated through a low-fidelity response model fl(x,t) as given by the blue surface. However, low-fidelity models typically include an additional set of model parameters t that need to be calibrated using a set of low-fidelity z and high-fidelity y data points [15]. These additional parameters are often introduced using approximation models (e.g., force field parameters associated with coarse-grained beads in molecular dynamics) or phenomenological models (e.g., constitutive law parameters in a finite element simulation). However, in many engineering challenges the model parameterst have a functional relation with the model variablesx (i.e., t=η(x) as represented by the yellow curve on the horizontal plane). For example, such challenges are encountered when calibrating the degradation of industrial equipment health indicators as a function of time [6,7], the calibration of interphase properties in composite materials as a function of the local structure [8], and the calibration of force field parameters in a coarse-grained polymer model with respect to temperature [9]. While a piece-wise linear approximation of the calibration functions can be achieved [9], this would also introduce a large number of parameters and typically results in unreasonably complex functions. Specifically, these piece-wise functions are likely to lead to erroneous results when evaluating the low-fidelity model with different initial conditions or length scales (e.g., a calibrated coarse-grained (CG) molecular dynamics simulation can be used to study physical phenomena at length and timescales that are too computationally exhaustive to be studied using the high-fidelity all-atomistic (AA) model). Consequently, we would therefore desire a statistical method to systematically balance the low-fidelity model complexity (i.e., reduce the number of model parameters) with the calibration accuracy (i.e., the low-fidelity model and high-fidelity model discrepancy). We address this challenge by proposing a structured semiparametric approach that helps modelers identify and calibrate an approximate set of low-dimensional deterministic calibration functions.

Fig. 1
Visualization of functional calibration as a nonisometric surface to curve matching. Through a set of high-fidelity data (light red circles on the right surface) and low-fidelity data (light blue circles in the plot area), we want to learn where the high-fidelity curve (red curve) lies on the low-fidelity response surface (blue surface) to learn the unknown calibration function (yellow line).
Fig. 1
Visualization of functional calibration as a nonisometric surface to curve matching. Through a set of high-fidelity data (light red circles on the right surface) and low-fidelity data (light blue circles in the plot area), we want to learn where the high-fidelity curve (red curve) lies on the low-fidelity response surface (blue surface) to learn the unknown calibration function (yellow line).
Close modal

In this paper, we make a distinction between model variables (x) that can be controlled by the modeler (e.g., system geometry, and material choice) and the model parameters (t) that are unique to the physical system being studied (e.g., Young's modulus, yield stress, and Poisson's ratio) [10]. In many cases, the true values of the model parameters are unknown as experiments are too expensive and/or the parameters cannot be measured directly. Calibration provides modelers with an approach to infer the values of these model parameters from two sets of data: (i) low-fidelity samples from a simulation model, and (ii) high-fidelity samples from either physical experiments or high-fidelity simulations [3,11]. This dataset can be used in a frequentist-based method that involves finding the values for the model parameters by minimizing the discrepancy between the low-fidelity and the high-fidelity models [11]. However, these methods generalize poorly as they do not provide confidence bounds on the identified model parameters. Alternatively, Bayesian approaches can address this limitation by quantifying the uncertainty in the model parameters [1214].

Parametric functional calibration can be achieved by assuming the class of calibration functions a priori (i.e., η(x;θ)) and then calibrating the unknown calibration parameters(θ) [15]. The experimental efficiency of functional calibration can be improved using emulators and the sequential acquisition of training samples [16]. In addition, the authors of Ref. [17] considered a variant of functional calibration where the low-fidelity and the high-fidelity models do not share the same set of model parameters. Specifically, they proposed a linear input mapping approach on the inputs of the low-fidelity model and used a regularization term to limit over-fitting. This allows for flexibility in the class of functions but cannot be used for functional calibration when the inputs to the high-fidelity model are also inputs to the low-fidelity model. If the functional forms of the model parameters are not known, then a prior distribution can be placed on the space of unknown calibration functions η, instead of the calibration parameters (θ), and updating them using the available training data [18,19]. This approach is typically referred to as nonparametric functional calibration. In addition, the calibration accuracy of nonparametric calibration can be improved by formulating and solving it as a “shortest path” problem [20,21]. Despite these contributions, the computational cost of functional calibration is often too exhaustive for problems with a large number of calibration functions, nor does it provide support for what class of functions should be used.

In this work, we introduce a semiparametric calibration method to select an appropriate class of calibration functions aimed at reducing the number of calibration parameters and limiting the compromise in prediction accuracy. We are specifically interested in the functional calibration of a low-fidelity model to an intrinsically noisy high-fidelity data source. Intrinsic noise is a type of aleatoric uncertainty that cannot be reduced with the acquisition of additional data points and is encountered in simulation models (e.g., agent-based models) as well as in physical experiments [22,23]. For example, such problems arise during the calibration of low-fidelity CG molecular dynamics simulations, to high-fidelity AA models. After finding the calibration functions, the modeler can use the CG model to explore physical phenomena at lengths and time scales that are intractable with the AA model. Note that changing the lengths and time scales is an invasive modeling procedure that precludes the modeler from using emulation techniques and emphasizes the importance of having an appropriate deterministic approximation of the calibration functions. An earlier work that addresses this issue involves the use of a global optimization method coupled with a neural network to minimize the discrepancy between the CG and AA model by tuning the force field calibration parameters [24]. However, a drawback of this approach is that the modeler needs to define the functional form of the force field parameters a priori, risking poor calibration accuracy and overfitting the calibration functions. The Bayesian-based method presented in this work circumvents this issue through the use of uncertainty quantification (UQ) to inform what calibration functions can be replaced with low-dimensional approximations.

The semiparametric calibration method introduced in this work provides support in selecting an appropriate class of functions for problems with a large number of calibration functions (i.e., η={η1,,ηnη}T where nη10). We consider this to be a semiparametric approach as a satisfying class of calibration functions is selected within the method as opposed to defining it a piori. This is similar to the decision-support framework presented in Ref. [25] to identify and remedy interpretability issues encountered when using Bayesian calibration. First, we assume a high-dimensional and flexible class of candidate functions for each model parameter. This introduces a large number of calibration parameters for which we quantify their uncertainty through a modified approximate Bayesian computation scheme. Specifically, we use the Hessian of a distance function at its global optimum to reduce the number of samples necessary to infer the first two statistical moments of the calibration parameters. The uncertainty sources considered in this method include the intrinsic uncertainty in the computer simulations as well as the interpolation uncertainty associated with the use of emulators. Through the posterior distribution of the unknown calibration functions, the modeler gains insight into the sensitivity of each calibration function, their signal-to-noise ratio, and their correlation. This insight is then used to identify what calibration functions can be replaced with low-dimensional approximations, and what class of functions should be used. Consequently, we present a systematic approach that facilitates the synthesis of multiple sources of domain knowledge to support modeling decisions (e.g., data from physics-based simulations, experimental data, and expert knowledge).

In the remaining sections of this paper, we first introduce available parametric and nonparametric functional calibration methods (Sec. 2). Next, we introduce the newly developed semiparametric functional calibration method and its constructs (Sec. 3). To verify the method's performance and validity, we continue our discussion by testing it on a test problem and an engineering example, with three and 14 unknown calibration functions, respectively (Sec. 4). Finally, we conclude this work by summarizing its findings and introducing opportunities for future work (Sec. 5).

2 Background

Model calibration enables synthesis of data from: (i) simulation models with controllable degrees of fidelity [2628], (ii) physical experiments and a single simulation model [12], and (iii) multiple simulation models [29]. In this section, we provide an introduction to the general formulation of available methods for model calibration and its extension to functional calibration.

2.1 Model Calibration.

Consider that we have a set of model variables xX that are the inputs to a low-fidelity simulation fl(·) and a high-fidelity simulation or experiment fh(·). Modelers are typically interested in learning how fl(·) and fh(·) depend on x, or to maximize fh(·) with respect to x. The data-efficiency of these objectives can be improved by the updating scheme proposed in Ref. [12] that relates the two data sources as
(1)

where n is the number of observed samples, η is a set of unknown but assumed to be constant modeling parameters, δ(·) is a bias function that captures the systematic difference between the low-fidelity and the high-fidelity data source, ρ is a scaling parameter, and ei is a noise term. The modeling parameters in this formulation are typically associated with some physical constants but cannot be observed experimentally. In addition, note that the high-fidelity function has been written as a function of the modeling variable x. Technically, it is also a function of the modeling constants η, but it is convention to not include this in the formulation as it cannot be observed experimentally, let alone controlled.

In Ref. [30], a calibration approach using the formulation presented in Eq. (1) is proposed that minimizes the distance between the low-fidelity and the high-fidelity model over a set of commonly observed variables xi,(i=1,,n) as
(2)

where ρ is assumed to be equal to 1, E is the space of admissible values of the modeling parameters, and ||·||L2 is the L2 norm. However, in the formulation of Eq. (2), the calibration parameters are assumed to be constant over the design space whereas many calibration problems show a functional relation with respect to the modeling variables [6,8,9].

2.2 Functional Model Calibration.

The available literature on functional calibration can be divided into two categories: parametric calibration, and nonparametric calibration. The parametric-based approaches assume the calibration functions to depend on the modeling variables as well as a set of calibration parameters θ (i.e., η(x,θ)). By defining the functional form a priori, the modeler can estimate the calibration parameters θ by minimizing the distance between the low-fidelity and the high-fidelity emulation models [31]. In contrast, an alternative line of study involves nonparametric approaches that substitute the modeling parameters η in Eq. (2) with a functional parameter η(x) (i.e., η=η(x)) [19]. As put forth in Ref. [16], we can replace the L2 norm with a sample average and add a regularization term to get a meaningful approximation of the calibration functions as
(3)
where nη is the number of calibration functions and λ is a smoothing parameter. In addition, the ηj is modeled on a reproducing kernel Hilbert space for which NCj is the native space generated from a kernel function Kj(·,·) with a corresponding norm ||ηj||NCj2 on the roughness of the jth calibration function [32]. By penalizing the roughness of the kernel functions, the algorithm places emphasis on smoother functions. In Eq. (3), we require that we query the low-fidelity and the high-fidelity data source an equal number of times at the same locations. This has the desirable property that the optimal solution lives in the finite-dimensional span of the training examples so that we can perform kernel operations independent from the infinite-dimensional feature space, as shown through the Generalized Representer theorem [33]. However, in most calibration scenarios the modeler can only query the low-fidelity and the high-fidelity data source a small number of times because of computational and experimental budgetary constraints. It is for this reason that the use of emulation models plays an important role in calibration as, they enable a designer to approximate the low-fidelity model fl(·) and bias function δ(·). A popular choice of emulators for calibration is Gaussian process models [28]. These types of emulators provide a posterior predictive distribution at any unobserved set of modeling variables, and calibration parameters for the low-fidelity data source, as
(4)

where N(·) is a normal distribution with mean μ(·) and standard deviation σ(·). In addition, y={y1,,yn}T is a set of observations associated with a training dataset of modeling variables (and calibration parameters for the low-fidelity data source). The GP is typically modeled through a squared exponential or a Matérn kernel function [34]. Consequently, GP emulators provide a practical approach that can be used to replace costly low-fidelity and high-fidelity simulations to expedite the calibration process.

3 A Semi-Parametric Method for Calibration

In this section, we introduce a method for calibration function formulation, UQ and dimension reduction.

3.1 Overview of the Semi-Parametric Calibration Method.

The flow of the semiparametric calibration method is shown in Fig. 2 and consists of three steps: (i) data acquisition and emulation, (ii) UQ of candidate functions, (iii) and reduced dimensional calibration. The calibration method starts with the creation of a design of experiments of uniformly distributed input conditions for the low-fidelity fl(·) and the high-fidelity model fh(·) (Step 1a). A design of experiments with nl samples for the low-fidelity model is given as {xl,tl}={{xl,1,tl,1}T,,{xl,nl,tl,nl}T}T where xl,iXd,(i=1,,nl) is a d-dimensional vector of model variables and tl,iTp,(i=1,,nl) is a p-dimensional vector of model parameters. Note that this means that the experimental design of the low-fidelity model takes samples from the joint space X×T. A design of experiments with nh samples for the high-fidelity model is given as xh={xh,1,,xh,nh}T where the vector xh,jXd,(j=1,,nh) is a d-dimensional vector of model variables defined on the same space as the low-fidelity model. Subsequently, we proceed by acquiring the response data for the low-fidelity model z=fl(xl,tl) and for the high-fidelity model y=fh(xh) by performing the necessary simulations. To predict the response for unobserved combinations of model variables xX, and modeling parameter tT we train two GP emulators as f̂l(x,t)|z and f̂h(x)|y for the low-fidelity model and high-fidelity model, respectively (Step 1b). Subsequently, we use these emulators in an optimization scheme (Step 1c) to find the deterministic optimal calibration parameters θ* for a candidate set with nη high-dimensional calibration functions η(x,θ)={η1(x,θ),,ηnη(x,θ)}T that we assumed to be radial basis functions (RBFs) due to their flexibility and ease of implementation [14,18,31].

Fig. 2
Flowchart for UQ-based decision support for semiparametric functional calibration. The method consists of three steps: (i)acquire training data to train emulation models, (ii) modified approximate Bayesian computation to quantify the uncertainty in the calibration parameters of the candidate set of high-dimensional calibration functions, and (iii) usage of the quantified uncertainty to simplify the calibration function parametrization.
Fig. 2
Flowchart for UQ-based decision support for semiparametric functional calibration. The method consists of three steps: (i)acquire training data to train emulation models, (ii) modified approximate Bayesian computation to quantify the uncertainty in the calibration parameters of the candidate set of high-dimensional calibration functions, and (iii) usage of the quantified uncertainty to simplify the calibration function parametrization.
Close modal

In contrast to finding a set of calibration parameters that minimizes the difference between the predicted mean responses, we want to infer the first two statistical moments of the calibration parameters as θ̂N(μ̂θ,Σ̂θ) using a distance metric D(θ|f̂h,f̂l,η). This lends itself to an approximate Bayesian computation inference scheme that involves sampling the posterior distribution of calibration parameters (Step 2a) without the need for evaluating the likelihood function. This is a significant advantage as evaluating the likelihood function requires the inversion of the covariance matrix that quickly becomes intractable for problems with a moderate to large number of calibration functions (nη10). In addition, we introduce a more data-efficient sampling process that uses the Hessian of the distance function to retrieve the eigenvectors of the calibration parameter covariance matrix (Step 2a and Step 2b). Subsequently, on each of these vectors, we find the two samples on the approximate Bayesian computation rejection boundary [35]. This enables a modeler to find two samples Xi,(i=1,,nθ) along each eigenvector that have the same distance function value D(·) (Step 2c). Subsequently, this provides a small yet efficient set of samples X={X1,,Xnθ}T from which we approximate the first two statistical moments of the candidate set of calibration parameters' posterior distribution as μ̂θ=E(X) and Σ̂θ=Cov(X), respectively (Step 2d). This assumes that the distribution of the calibration parameter is jointly Gaussian, and as such, we used a hat to indicate that this is an approximation.

The posterior distribution of the calibration parameters θ̂N(μ̂θ,Σ̂θ) can be used to draw insight into the functional form of the calibration functions. We achieve this by visualizing the posterior distribution of the calibration functions η(·) and by evaluating their signal-to-noise ratio (Step 3a). Subsequently, the confidence bounds, as well as the covariance matrix, provide insight into the underlying class of functions (Step 3b). If the confidence bounds are narrow and the mean function shows a strong trend (i.e., the expected range of the function is large) then this informs the modeler that the calibration functional forms need to be selected carefully. Moreover, in this scenario, the signal-to-noise ratio will be relatively large and a more complex calibration functional form (i.e., more input dimensions) can be selected. Conversely, if the confidence bounds are wide, or the trend is weak, then the signal-to-noise ratio will be small and the functional form of this calibration function is less consequential to the calibration accuracy and can thus be simplified (e.g., constant, or linear) without significantly affecting the calibration accuracy. Finally, the simplified set of calibration functions η(x,ϑ)={η1(x,ϑ),,ηnη(x,ϑ)}T can then be used to find the optimal set of reduced dimensional calibration parameters ϑ*. The presented method helps a modeler reduce the complexity of the calibration functions, provides insight into the physics that underlie their models, provide strong calibration accuracy, and enables the generalizability of these calibration functions to low-fidelity simulations at larger length and time scales. In the remainder of this section, we will provide additional details on the method's individual constructs: (i) distance and calibration functions, (ii) modified approximate Bayesian calibration scheme, and (iii) how to select low-dimensional calibration functions.

3.2 Formulation of Distance and Calibration Functions.

We draw inspiration from recent work on the use of the stochastic gradient descent optimization trajectory to approximate the first two statistical moments of a neural network's model parameters [36,37]. The assumption in Ref. [37] is that the curvature of the loss function around the global optimum provides meaningful information about the covariance matrix of a neural network's parameters. In our case, we need to define a loss function that adequately captures the goodness of fit that includes the mean prediction as well as the interpolation uncertainty. For this reason, we use GP emulators as they provide an approximation of the posterior distribution at unobserved spatial locations. Specifically, the emulator as defined in Eq. (4) provides a posterior predictive distribution that gives the probability of observing response f from our low-fidelity model and high-fidelity model as P(f̂l(x,η(x,θ))=f|z) and P(f̂h(x)=f|y), respectively. Observe that both the low-fidelity and the high-fidelity emulators predict the same response, and so we can evaluate the goodness of fit for a specific set of calibration parameters θ and calibration functions η at a specific model variable x˜ as
(5)

where F is the range of the response variable. Note that in this expression we calculate the area under the curve of the probability predicted from the product of the low-fidelity emulator with the high-fidelity emulator. This has two desirable properties: (i) minimizing the distance function will force the mean predictions of the emulators to be the same, and (ii) the distance function places a higher importance on improving the goodness of fit in areas of the response surface with lower posterior predictive uncertainty. Another interesting observation is that the range of Eq. (5) is D(θ)(0,1], where D(θ)=0 would correspond to a θ such that the two GP emulators are equivalent (i.e., μh(x˜)=μl(x˜,η(x˜,θ)) and σh(x˜)=σl(x˜,η(x˜,θ))=0). However, in most cases the data observed from the high-fidelity model is corrupted with noise (i.e., σ̂h(x˜)>0) and as such the lower bound will be greater than zero. This can be corrected by taking the root of the term inside the integration of Eq. (5) and would provide a probabilistic interpretation for the distance function. However, doing so does not alter the calibration performance and thus we leave it as an interesting observation.

From Eq. (5), we can approximate the goodness of fit for of a specific set of calibration parameters θ at a fixed value of the model variables. Consequently, by integrating over the range of the model variables we obtain the following distance function
(6)

If we normalize the range of the model variables (i.e., Θ[0,1]d), then the aforementioned probabilistic interpretation of the distance function is retained.

With the distance function defined, we need to adopt a flexible class of functions that can be used as a candidate set of calibration functions η(x,θ). For this purpose, we suggest the use of RBFs as they provide a highly flexible class of functions that are easy to implement. The GP emulators introduced in Sec. 2 are a special case of RBF functions in that their kernel is symmetric and positive definite. The general formulation of an RBF function with nc centers and a squared exponential kernel as
(7)

where the ith element of the nc×1-dimensional vector k is defined as exp(wT(xi(c)x)2), the i, jth element of the nc×nc matrix K is defined as exp(wT(xi(c)xj(c))2), and c is an nc×1 vector of weights (i.e., a set of parameters that need to be inferred from the training data). In this expression, the terms kTK1 provide a set of nc radial functions centered at spatial locations x(c). Taking the inner product of these radial functions with the weights c allows us to model a wide range of continuous functions. The nc center locations x(c)={x1(c),,xnc(c)}T are best designed through a space-filling set of centers, in this work, we will use a uniform grid (i.e., a grid with a total of ncnηi centers, with nc centers for each of the nηi model variables that are an input to the ith calibration function). This introduces a set of calibration parameters θ={w1T,,wnηT,c1T,,cnηT}T to parameterize a total of nη calibration functions. The form of the calibration functions is equivalent to those used in Ref. [38], except that we use a smaller number of center locations. The big advantage of using RBF functions is that we can readily control the complexity of the calibration functions that we can capture by increasing the number of centers. In addition, using the RBF functions we can easily avoid extrapolating our emulation model f̂l(·) by setting the range for the centers equivalent to the range of the associated dimension in the training dataset. In contrast, setting the right range for the constants of a set of polynomials is less straightforward [39].

3.3 Uncertainty Quantification of the Calibration Functions.

Using the distance metric in Eq. (6), we now want to directly sample the posterior distribution of the calibration functions through an approximate Bayesian computation scheme. This can be achieved through available methods but is inefficient for high-dimensional problems [40,41]. Consequently, we introduce an alternative approach that requires fewer samples of the posterior distribution of the calibration parameters. In Eq. (6), the GP emulators for the low-fidelity model (f̂l(·)) and the high-fidelity model (f̂f(·)) are used to create a set of pseudo-observations for any value of the calibration parameters θ [35]. Two typical approaches to sampling the posterior distribution of the calibration functions is the approximate Bayesian computation rejection algorithm [40] and the approximate Bayesian computation Markov chain Monte Carlo algorithm [41]. However, these approaches are known to be inefficient and intractable for problems with a large number of parameters [35,42]. In this subsection, we introduce a modified version of the approximate Bayesian computation rejection algorithm that reduces the number of necessary samples to approximate the first two statistical moments of the calibration function's posterior distribution. Specifically, the algorithm starts by identifying the eigenvectors of the posterior distribution of the calibration parameters and then finds the two samples over each vector that are at the rejection boundary.

The approximate Bayesian computation rejection scheme involves sampling from the prior distribution of calibration parameters P(θ) and then retaining these samples if D(θ|f̂h,f̂l,η)<ξ, where ξ is a rejection criteria identified by the modeler. The outcome of this procedure is a set of samples that are approximately distributed as the posterior distribution of the calibration parameters. The first step in the calibration method is to obtain a deterministic approximation of the optimal calibration parameters (Step 1c) as
(8)

The optimization of this function can be achieved through any gradient-based method starting from multiple initial calibration parameters or an evolutionary-based algorithm [43].

Using the optimal set of calibration parameters θ* we can approximate the eigenvectors of the posterior distribution from the Hessian matrix H. Specifically, the i, jth element of the Hessian matrix is obtained as
(9)

Subsequently, the eigenvectors can be obtained through an eigenvalue decomposition as H=QΛQ1, where the nθ columns of Q are the eigenvectors qi,(i=1,,nθ), and Λ is a matrix of mostly zeros except for its diagonals that contain the eigenvalues. Using the Hessian to infer the eigenvector is similar to taking the quadratic term of the local Taylor's series expansion at θ* [43]. Moreover, it should be noted that using a Newton-type gradient-based optimization scheme (e.g., sequential quadratic programming) to find the global minimum of the distance function has the added advantage of providing an approximation of the Hessian H.

From the set of eigenvectors, we now need to infer the posterior distribution of the calibration functions (Step 2c and 2d). This can be done by finding the two samples with the same distance function values over each eigenvector. This idea has been visualized in Fig. 3 where the red lines indicate the eigenvectors that have been identified in Steps 2a and 2b and the surface plot represents the distance function given in Eq. (6). The algorithm then searches for two points Xi,(i=1,,2nθ) over each eigenvector that have a distance function value equal to ξ=D(θ*)/2. In other words, we have selected the rejection boundary of the value to be half of the distance function response at the global optimum μ̂θ. This fraction can be changed but will mainly influence the magnitude of the covariance matrix. This is because the equal distance function boundary approximates an ellipsoid (as can be seen from Fig. 3). Note that each element of the sample set X={X1,,Xnθ}T can be found by solving a one-dimensional optimization problem. Subsequently, we can approximate the first two statistical moments of the posterior distribution of the calibration parameters as μ̂θ=E(X) and Σ̂θ=Cov(X).

Fig. 3
Identification of samples at the rejection boundary of the calibration parameter covariance matrix. The equal distance function samples (green circles) are obtained by performing a one-dimensional optimization over each eigenvector (red lines).
Fig. 3
Identification of samples at the rejection boundary of the calibration parameter covariance matrix. The equal distance function samples (green circles) are obtained by performing a one-dimensional optimization over each eigenvector (red lines).
Close modal

Note that this approach is a more data-efficient alternative to conventional approximate Bayesian computational rejection algorithm that would require significantly more samples from the prior distribution to approximate the posterior distribution. More specifically, by drawing samples from a noninformative prior distribution the conventional rejection algorithm would find a set of samples uniformly covering the elliptically shaped surface of the distance function in Fig. 3. Consequently, from this observation, we can conclude that the introduced sampling scheme will find a covariance matrix that would be the same as the approximate Bayesian computation rejection scheme with a uniform/uninformative prior. This implies that the introduced scheme does not enable a modeler to use an informative prior distribution for the calibration parameter, but we believe that this type of information will only rarely be available.

An analogy to the procedures in Steps 2a to 2d is that we are trying to find the minimum distance function value for which we can fit an nθ-dimensional hyper-ellipse with a constant surface area over the distance function surface D(·). This is an optimization problem that is remarkably similar to robust design optimization where modelers are concerned with finding a design that is robust with respect to various sources of uncertainty [44,45]. However, the difference with our approach is that the surface area of this hyper-ellipsoid is implicitly selected through the rejection boundary a priori.

3.4 Approximating the Class of Calibration Functions.

In the final steps of the method (Step 3a to 3c), we use the posterior predictive distribution of the unknown calibration functions to inform their parametric form. This process starts with propagating the uncertainty in the calibration parameters θ̂N(μ̂θ,Σ̂θ) to quantify the uncertainty in the calibration functions η(·). The quantified uncertainty and the range of the function can then be used to inform their parametrization. For example, consider the two plots in the bottom left of Fig. 2, here we have two functions with roughly the same magnitude of confidence bounds, but one shows a strong linear upward trend whereas the other appears to be constant. Consequently, the modeler has reason to believe that these functions can be approximated to be linear and constant, respectively. The simplified set of calibration functions η(·) can be established by considering multiple classes of low-dimensional functions (e.g., constant, linear, quadratic, sigmoid, and RBF). By comparing the relative difference between the uncertainty bounds and the shape of the candidate set of calibration functions the modeler can observe the most likely shape of each function and their relative importance to the calibration accuracy. However, visual inspection is not systematic and becomes complicated for calibration functions with more than two inputs (i.e., nηi>2).

When assessing the posterior distribution of the calibration functions we can make two observations: (i) the range of the unknown functions (i.e., signal), and (ii) the magnitude of the uncertainty (i.e., noise). Consequently, it would be interesting to look at the signal-to-noise ratio Si,(i=1,,nη) that we define as
(10)

where the numerator is the expected range of the ith unknown calibration function and the denominator is the maximum posterior standard deviation. A low signal-to-noise ratio (e.g., S3) indicates poor identifiability as the range of the function is small compared to the range of the confidence bounds. Consequently, these functions can be simplified significantly as their form is of marginal consequence to the model's prediction accuracy. In contrast, if the signal-to-noise ratio is large (e.g., 6<S), then this implies that the training data provides significant insight into the unknown calibration functions. Consequently, it is reasonable to select a more complex functional form that provides a good approximation of the unknown function. Finally, for calibration functions with a moderate signal to noise ratio (e.g., 3<S6) it is possible to glean some insight into their functional form. However, these functions are best described through low-dimensional approximations (e.g., linear or quadratic) as the “true” functional form is still likely obscured by the nontrivial confidence bounds. The distinction of the small, medium, and large signal-to-noise ratios has been established based on empirical results as presented in Sec. 4. A summary of recommendations to guide the selection of the appropriate functional forms is provided in Table 1. Alternatively, variance-based sensitivity analysis can be used to arrive at a metric to inform the class of calibration functions [46]. The advantage of variance-based sensitivity analysis is that it provides an approximation of the higher-order interaction effects. However, in this work we recommend the use of the signal-to-noise ratio as variance-based sensitivity analysis can be time intensive for high-dimensional problems and the interaction affects can also be gleaned from the off-diagonal element of the calibration parameter covariance matrix Σ̂θ.

Table 1

Summary of recommendations for the selection of an appropriate class of calibration functions with respect to the inferred signal to noise ratio

Signal to noise ratio SSmall signal to noise 0<S3Medium signal to noise 3<S6Large signal to noise 6<S<
Calibration functional form recommendationGreat potential for simplification, a constant function will be appropriate in most cases.Some potential for simplification, a linear or quadratic function will be appropriate in most cases.The calibration function will have to be selected carefully to match the semiparametric approximation.
Additional considerationsThis scenario indicates that the calibration functions are unidentifiable if the noise is relatively large.These functions can be simplified with low-dimensional approximations and provide some insight into their “true” form.The general trend of the approximate calibration function is likely consistent with the “true” form.
Signal to noise ratio SSmall signal to noise 0<S3Medium signal to noise 3<S6Large signal to noise 6<S<
Calibration functional form recommendationGreat potential for simplification, a constant function will be appropriate in most cases.Some potential for simplification, a linear or quadratic function will be appropriate in most cases.The calibration function will have to be selected carefully to match the semiparametric approximation.
Additional considerationsThis scenario indicates that the calibration functions are unidentifiable if the noise is relatively large.These functions can be simplified with low-dimensional approximations and provide some insight into their “true” form.The general trend of the approximate calibration function is likely consistent with the “true” form.

4 Results and Observations

The use of the introduced calibration method will be highlighted in this section by applying it to a test problem with three calibration functions, and the calibration of a low-fidelity CG simulation to a high-fidelity AA simulation with 14 unknown calibration functions.

4.1 Four-Dimensional Test Problem.

To demonstrate the functionality of the proposed method, we will start by applying it to a test problem that has one modeling variable and three calibration functions. Specifically, the high-fidelity model is defined as
(11)
The low-fidelity function is similar to the high-fidelity function except that the calibration functions η(·)={η1(·),η2(·),η3(·)}T are unknown, and we run the simulation considering them as scalar inputs
(12)

where t={t1,t2,t3}[1.5,1.5]×[0.5,1.5]×[0.25,0.75]3. In order to make inferences into the form of the assumed to be unknown calibration functions η, we start by creating a training dataset (Step 1a) including 100 samples of the high-fidelity model and 400 samples of the low-fidelity model. The samples of the low-fidelity model span a four-dimensional space that consists of the model variable x as well as the calibration parameters t. Subsequently, we proceed by training our emulation models f̂h(x)|yN(μ̂h(x),σ̂h(x)) and f̂l(x,t)|zN(μ̂l(x,t),σ̂l(x,t)) (Step 1b) and then use these to find the deterministic optimal set of parameters θ* of the candidate set of calibration functions (i.e., using Eq. (6) as part of Step 1c). Moreover, we have set the number of RBF centers equal to five (i.e., nc = 5).

Using the approximate Bayesian computation scheme (Step 2) we acquire the posterior distribution of the calibration parameters θ̂N(μ̂θ,Σ̂θ) and use it to get the posterior distribution of the calibration functions. Specifically, by drawing 1e4 samples, we retrieve the calibration functions as presented by the purple shaded regions and the dashed-dotted lines in Figs. 4(a)4(c). From these figures it can be observed that η1(·) and η2(·) have relatively small confidence bounds compared to η3(·). Moreover, this gives us the signal-to-noise ratios as S1=14.01,S2=3.21, and S2=1.25. This informs us that selecting an appropriate functional form for η1(·) and η2(·) is more important than for η3(·). This observation is consistent with our expectation, as when we look at the high-fidelity function (Eq. (12)) we find that the third calibration function influences the frequency of a sinusoidal term and, as a result, variations in η3(·) will decrease the value of Eq. (6) only marginally. More specifically, by observing the predicted response of the high-fidelity model and the low-fidelity model in Figs. 4(e) and 4(d), we find that the influence of the sinusoidal trend associated with the third calibration function is marginal. Moreover, if the frequency of this sinusoidal trend were to change as a result of η3(·), then, the confidence bounds will still overlap. However, changes in η1(·) or η2(·) are proportional to the deviation between the low-fidelity model and the high-fidelity model.

Fig. 4
Calibration results of a one-dimensional function with three unknown calibration functions. The uncertainty in the calibration function has been approximated using three RBFs with five centers and enables the simplification of their functional form with a limited compromise to the calibration accuracy: (a) Quantified uncertainty for the first calibration function, the simplified functional approximation, and the true function. (b) Quantified uncertainty for the second calibration function, the simplified functional approximation, and the true function. (c) Quantified uncertainty for the third calibration function, the simplified functional approximation, and the true function. (d) Comparison of the high-fidelity with the low-fidelity response for the most probable calibration parameters. (e) Comparison of the high-fidelity with the low-fidelity response for the reduced dimensional calibration functions.
Fig. 4
Calibration results of a one-dimensional function with three unknown calibration functions. The uncertainty in the calibration function has been approximated using three RBFs with five centers and enables the simplification of their functional form with a limited compromise to the calibration accuracy: (a) Quantified uncertainty for the first calibration function, the simplified functional approximation, and the true function. (b) Quantified uncertainty for the second calibration function, the simplified functional approximation, and the true function. (c) Quantified uncertainty for the third calibration function, the simplified functional approximation, and the true function. (d) Comparison of the high-fidelity with the low-fidelity response for the most probable calibration parameters. (e) Comparison of the high-fidelity with the low-fidelity response for the reduced dimensional calibration functions.
Close modal

Comparing Fig. 4(a) with Fig. 4(b) we observe that the mean trend of the first calibration function (i.e., μη1(·)) changes more significantly than the mean trend of the second calibration function (i.e., μη2(·)). Subsequently, we have evidence that it is more reasonable to assume a more high-dimensional and complex function for the first calibration function (i.e., μη1(·)). This is an encouraging observation as this is consistent with the relatively high signal-to-noise ratio and insights that we can draw from the problem formulation as stated in Eq. (12) (i.e., we know that the “true” functional form is logistic). Modeling such a function rigorously would require at least four calibration parameters. In contrast, η̂2(·)μη2 shows a much simpler linear trend. Nevertheless, the confidence bounds can hide more complex functions, the only conclusion that we can make from this observation is that we can model the second calibration function with a linear class of functions and still expect reasonable calibration accuracy. This observation is even more relevant for the third calibration parameter where the confidence interval could hide a simple constant function or a wildly fluctuating multimodal function that would require many variables to parametrize. Nevertheless, these figures inform us that we can expect to have a reasonable calibration accuracy if we simplify the functional class of η1(·),η2(·), and η3(·) by an RBF with five centers a linear function and constant function, respectively. Moreover, note that these conclusions are consistent with the recommendations presented in Table 1.

With the simplified functional form of the three calibration functions ηi,(i=1,2,3), we can once again minimize Eq. (6) to get the reduced-dimensional set of optimal calibration parameters ϑ*. This is a significant simplification in the number of parameters, as the RBF calibration with five centers for three calibration functions requires 18 parameters whereas the simplified set of functions only requires only nine parameters. In addition, the method allows a modeler to make an informed decision in terms of the class of functions to assume for each calibration function. If a modeler was to consider using a set of five functional classes (e.g., constant, linear, quadratic, cubic, RBF), then a problem with three calibration functions would have a total of 53=125 alternative formulations. Exploring all these alternatives manually would be a tedious and time-consuming effort.

A concern that arises in the three-dimensional test problem is that the simplification will have a negative influence on calibration accuracy. When plotting the predicted response in Fig. 4(d) of the low-fidelity model and the high-fidelity model, we observe that it is difficult to separate the two from one another when using the high-dimensional functional parameterization (i.e., the dashed-dotted purple lines in Figs. 4(a)4(c)). Contrasting this with that of the simplified functionalization (i.e., the solid black lines in Figs. 4(a)4(c)) in Fig. 4(e), we find that the visual difference is hard to notice. In fact, the root mean squared percentage error (RMSPE) [47,48] for the original functionalization and the simplified functionalization are 0.36% and 0.38%, respectively. As might be expected, the calibration accuracy of the simplified functionalization is slightly worse, but still very good. In fact, an RMSPE less than 10% is typically regarded as excellent [49]. This is further corroborated by the observation that it is difficult to differentiate between Figs. 4(d) and 4(e).

In order to compare the tradeoff between the reduction in model complexity and calibration accuracy we can compare the results using the Bayesian information criterion that is given as BIC=nlog(nθ)2log(L), where L is the likelihood and nθ is the number of calibration parameters [50]. By calculating the Bayesian information criterion for observing the high-fidelity training dataset by sampling from the low-fidelity emulation model, we get BIC =673.05 and BIC=615.15 for the mean approximation of the candidate set of functions and the simplified functions, respectively. The significance of this improvement should not be underestimated as it is on a logarithmic scale. Even though the Bayesian information criteria provide a statistically rigorous tradeoff between calibration accuracy and functional complexity, using it directly to optimize the objective would be complicated. This is because we would have to search over a design space with nfnη categorical alternatives, where nf is the number of considered functional classes.

4.2 Calibration of a Coarse-Grained Epoxy Model.

Modelers that are interested in exploring the thermomechanical properties of new epoxy-based composite materials frequently rely on simulation models to explore the high-dimensional design space of atomic structures and chemistries. The simulation models of choice are typically molecular dynamics; however, because AA simulations are intractable at experimentally verifiable length scales, modelers rely on CG simulations to predict thermomechanical properties. The force field of AA simulations is often well-defined, but CG models use empirical force fields with interaction parameters that must first be calibrated, for example, by matching the predictions of the higher-fidelity AA simulations. Moreover, previous studies have shown that the degree of crosslinking (DC) has a consequential effect on the properties of an epoxy system [5156]. Despite this knowledge, previous efforts have been unable to discover the functional dependence of the model parameters on the DC. In this study, we will use the introduced calibration method to elucidate this functional dependence. The material science related findings of this study have already been published in a previous publication [57]. In contrast, in this paper, we demonstrate the design related findings.

We use an epoxy compound with Bisphenol A diglycidyl ether (DGEBA) as the resin and either 4,4-Diaminodicyclohexylmethane (PACM) or polyoxypropylene diamines (Jeffamine D400) as the curing agent. AA simulations are used as our high-fidelity model and CG simulations are used as our low-fidelity model. Moreover, we consider DC as the modeling variable, and each unique coarse-grained bead has two calibration functions, cohesive energy, and bead size. In the coarse grain procedure, we identified seven unique beads for both the PACM and D400 systems, and thus, we have a total of 14 calibration functions. We consider these calibration functions to be the inputs to the low-fidelity CG simulation that enables us to predict the modulus, yield strength, Debye-Waller factor (u2), and the density of the two material systems (DGEBA + PACM and the DGEBA + D400). These properties have been selected as they are important markers of the overall thermomechanical properties. Using the introduced calibration method, we will simultaneously match all eight low-fidelity simulation model responses to the AA model responses. To calibrate the unknown functions with respect to multiple responses simultaneously we take their product at the inner integral of Eq. (6).

We started the calibration procedure by running 700 and 500 low-fidelity simulations for the DGEBA + PACM and DGEBA + D400 system, respectively. These simulations were performed using a fully sequential space-filling design obtain from a Sobol sequence [58]. The advantage of using a fully sequential space-filling design over a one-shot space-filling design is that we can sequentially add more simulation experiments until an emulation model with sufficient fidelity has been obtained [59,60]. In fact, this is the reason why 700 simulations where used to approximate the DGEBA + PACM CG model as opposed to the 500 used for the DGEBA + D400 system. In addition, we used a similar procedure for the high-fidelity model that included 19 simulations and 20 simulations for the DGEBA + PACM and DGEBA + D400 system, respectively. To quantify the uncertainty in the calibration parameters (Steps 2a to 2d) we assumed that for each of the 14 RBF functions we used three center points and one shared length-scale parameter. Subsequently, we have a total of nθ=43 calibration parameters for which we approximated their first two statistical moments. Using the approximate distribution of the calibration parameters, we quantified the uncertainty in the calibration functions as plotted in Fig. 5. In these plots, the purple shaded regions are the 95% confidence intervals and the purple dashed-dotted lines are the mean calibration functions. Moreover, the interaction energies are associated with the first seven calibration functions (i.e., εi=ηi,(i=1,,7)) and the bead sizes are associated with the last seven calibration functions (i.e., σi=η7+i,(i=1,,7)). From these figures, we can infer that many of the calibration functions can be replaced with either a constant or linear functions. Specifically, using the signal-to-noise ration for each calibration function (S1=1.77,S2=3.22,S3=6.22,S4=0.95,S5=3.54,S6=1.05,S7=0.84,S8=4.54,S9=2.87,S10=3.01,S11=1.83,S12=0.042,S13=1.62,S14=2.36) we identified selected the functional form according to the followings conditions, if 0<S3 we used a constant function, for 3<S6 we have used a linear approximation, and for 6<S we have used an RBF with three centers. The selected and calibrated functions, as inferred from the training data, are shown by the black solid lines in Fig. 5. Finally, some simplified calibration functions display a nontrivial difference between the quantified confidence bounds. While this observation highlights the high-dimensional/complex interaction between the calibration functions, it also confirms that the simplified functions have a form that is consistent with the form of the all-RBF approximations.

Fig. 5
Calibration functions quantified uncertainty and simplified approximation for a CG simulation. Purple dash-dotted lines represent the high-dimensional all-RBF approximation of the calibration functions, with the quantified uncertainty represented through the purple-shaded region. In addition, the solid black lines represent the simplified parametrizations informed from the quantified uncertainty.
Fig. 5
Calibration functions quantified uncertainty and simplified approximation for a CG simulation. Purple dash-dotted lines represent the high-dimensional all-RBF approximation of the calibration functions, with the quantified uncertainty represented through the purple-shaded region. In addition, the solid black lines represent the simplified parametrizations informed from the quantified uncertainty.
Close modal

The mean and the confidence interval for the eight CG model responses as predicted from the GP emulators for the simplified set of calibration parameters ϑ* are plotted by the blue shaded region and lines in Fig. 6, respectively. The orange plots represent the high-fidelity simulation response. From these figures, we can observe that the response of the low-fidelity simulation model matches well with that of the high-fidelity response. In particular, we observe that the mean trend of the density appears to provide the best match to the AA model. The reason for this is that in the distance function in Eq. (6), we consider the predictive uncertainty of the GP emulators and as such the algorithm places more emphasis on those responses with smaller intrinsic modeling uncertainties.

Fig. 6
Response functions for the CG model and the AA model given the all-RBF approximation of the calibration functions. The uncertainty bounds of the CG model have inflated as a consequence of using the distribution of the high-dimensional all-RBF calibration functions.
Fig. 6
Response functions for the CG model and the AA model given the all-RBF approximation of the calibration functions. The uncertainty bounds of the CG model have inflated as a consequence of using the distribution of the high-dimensional all-RBF calibration functions.
Close modal

In terms of fidelity, it would be of interest to compare the calibration accuracy of the simplified parametrization ϑ* with the mean approximation of the all-RBF parametrization μ̂θ, and thus we have plotted the responses in Fig. 7. From these plots, we observe that the candidate set of all-RBF parametrization performs significantly better in terms of Young's modulus and the Yield stress. However, this comes at the expense of accuracy in the predicted Debye-Waller factor, (u2). When comparing the accuracy of the mean prediction we find that the simplified parametrization and the candidate set of calibration functions have an RMSPE of 10.0% and 9.8%, respectively. As expected, the simplified parameterization performs slightly worse, but this change is only marginal compared to the decrease in the number of parameters necessary to parametrize the calibration functions.

Fig. 7
Response functions for the CG model and the AA model given the deterministic approximation of the reduced dimensional set of calibration functions. The mean response of the density functions has a closer match as a consequence of the relatively lower intrinsic modeling uncertainty.
Fig. 7
Response functions for the CG model and the AA model given the deterministic approximation of the reduced dimensional set of calibration functions. The mean response of the density functions has a closer match as a consequence of the relatively lower intrinsic modeling uncertainty.
Close modal

The simplified set of calibration functions contains 1 RBF, 4 linear, and 9 constant functions for a total of 21 calibration parameters. This is a significant reduction compared to the 43 parameters from which we started. In addition, if we were to consider that each calibration function could be modeled through one of five functions then the modeler would have to choose one class of functions from a set of 514 different calibration function combinations. Exploring more than 6 billion alternative sets of calibration functions is too time-consuming and as such the introduced calibration method provides a modeler with valuable support. The simplified set of calibration functions have been selected to reduce the number of calibration parameters and avoid over/under-fitting so that they can be used to reliable explore the properties of larger material systems using the low-fidelity CG model.

A final insight that can be obtained from this analysis is that the average absolute correlation between the calibration parameters of different calibration functions as plotted in Fig. 8. The conclusion that can be drawn from this figure is that η3(·) and η5(·) have a high correlation with many other calibration functions. This observation is in line with what is expected as these calibration functions are associated with the cohesive energy of the coarse-grained beads that are found in both the DGEBA + PACM and the DGEBA + D400 system. In addition, we can observe that the correlation between the calibration functions associated with the size of the beads (i.e., η8,,η14) are correlated among themselves but have little correlation with the interaction energies. This should come as little surprise as the density mainly depends on the bead size and the other mechanical properties depend on the cohesive energy. We make this observation here as it validates the potential of the introduced calibration method, for a more detailed discussion on the modeling aspects of this study we refer the reader to a paired study [57].

Fig. 8
The average absolute covariance between the RBF centers of the calibration functions for a coarse-grained epoxy model. The correlations illustrate that function three and function five have a high correlation with many of the other calibration functions. This is an encouraging observation as these calibration functions are associated with the coarse-grained beads that are found in both the DGEBA + PACM material system and the DGEBA + D400 material system.
Fig. 8
The average absolute covariance between the RBF centers of the calibration functions for a coarse-grained epoxy model. The correlations illustrate that function three and function five have a high correlation with many of the other calibration functions. This is an encouraging observation as these calibration functions are associated with the coarse-grained beads that are found in both the DGEBA + PACM material system and the DGEBA + D400 material system.
Close modal

4.3 Validation of Selected Calibration Functional Form.

A limitation of the introduced semiparametric calibration scheme is that it requires the modeler to select the appropriate calibration functional form based on the signal-to-noise ratios. As a consequence, given the same training dataset, different modelers can arrive at different sets of calibration functions η(·). We argue that this is an advantage of the proposed method, as it enables the modeler to leverage their implicit knowledge when selecting the calibration functional form. However, a counterargument would be that in the ideal scenario, this knowledge is explicitly accounted for by the physics of the simulation models and that the calibration functions are selected autonomously. Consequently, to provide additional corroboration of the introduced semiparametric calibration method, it would be beneficial to compare it to an approach where the selection of the calibration functional form is automated.

The autonomous selection of the calibration functional form is a challenging problem that requires the search over a combinatorial space of alternative functions. One approach to simplifying this problem is by considering only polynomial functions for the calibration functions and regularizing higher-order terms. Doing so considering only constant and linear terms of a polynomial we get the following calibration problem
(13)

where the calibration functions are given as η(xj,β)=β0+β1xj=[β0,1+β1,1xj,,β0,nη+β1,nηxj]T and yi,j is the ith response of the jth high-fidelity observations. Note that in Eq. (13) the regularization is only with respect to the parameters associated with the linear terms of the polynomials (i.e., β1=[β1,1,,β1,nη]T). Consequently, by minimizing Eq. (13) we only penalize functions that deviate from being constant. The parameters are regularized with respect to the L2 norm, the advantage of this is that we have a continuous objective function. However, the disadvantage is that in the optimization scheme parameters that are close to zero have a small gradient, and thus might not converge to zero as they would using the discontinuous L1 norm. Finally, we believe that this is an apt approach for validation as minimizing Eq. (13) can be expected to be on par with the introduced semiparametric calibration approach. The reason for this is that in the final set of calibration functions only one was not linear or constant. While the degree of model simplification is problem dependent, the primary advantage of the proposed method is that it informs when simplification of the calibration functions is realistic (i.e., when the signal-to-noise ratio is small).

A feature of Eq. (13) is that we can control the tradeoff between calibration accuracy and model sparsity through the regularization term γ. Setting a small value for γ will minimize the penalty on nonconstant calibration functions and thus provide more flexibility. In contrast, selecting a large value for γ will emphasize model sparsity by penalizing nonconstant calibration functions. To investigate the influence of the regularization term, we have calibrated the CG model using Eq. (13) for uniformly spaced values of γ{0,,5}. The resulting calibration accuracy has been plotted in Fig. 9(a) and the model complexity has been given in Fig. 9(b). The model complexity is the sum of the number of constant terms and the number of linear terms with a value greater than 1% of the admissible range of the calibration functions. Specifically, 1|β1,in|0.01 is an identity function where β1,in is the ith normalized nonlinear parameter. The discrete calibration results for 100 values of γ are presented by the blue circles and the blue lines present a moving average with a moving window of five. For comparison we also included the results obtained from the introduced semiparametric calibration approach, see the red dashed line. In addition, the black dotted line represents the first set of calibration functions obtained from Eq. (13) that was equal in model complexity to the results obtained from the introduced semiparametric calibration approach.

Fig. 9
Comparison of model complexity versus calibration accuracy through regularization. Model complexity ranges from 14 parameters (all constants) to 28 parameters (all linear functions). For comparison, we have included the results of the UQ informed calibration function selection by the red dashed lines. In addition, the autonomously selected calibration functions with equal complexity (14 parameters) are given by the dotted black lines. (a) Calibration accuracy (RMSPE) for different regularization weights γ. (b) Model complexity (number of parameters) for different regularization weights γ.
Fig. 9
Comparison of model complexity versus calibration accuracy through regularization. Model complexity ranges from 14 parameters (all constants) to 28 parameters (all linear functions). For comparison, we have included the results of the UQ informed calibration function selection by the red dashed lines. In addition, the autonomously selected calibration functions with equal complexity (14 parameters) are given by the dotted black lines. (a) Calibration accuracy (RMSPE) for different regularization weights γ. (b) Model complexity (number of parameters) for different regularization weights γ.
Close modal

It is encouraging to observe that only considering constant and nonlinear calibration functions nearly doubles the RMSPE compared to the introduced semiparametric calibration method. Moreover, it shows that considering more complex functions is a valuable exercise as it can greatly improve calibration accuracy without compromising model complexity. Specifically, from these observations we can make two conclusions:

  1. Decreasing γ increases model complexity and results in improved calibration accuracy. However, the calibration accuracy is not reduced to zero and as such considering only constant and linear functions does not provide sufficient flexibility.

  2. Increasing γ decreases model complexity and calibration accuracy. However, the model complexity does not reduce to zero, and as such it would be valuable to use an L1 norm in Eq. (13) to truly enable the selection of the functional form.

Nevertheless, the introduced semiparametric calibration method provides a valuable method to balance model complexity with calibration accuracy albeit dependent on the modeler's intuition or expert knowledge.

From the previous study, we can also learn the importance of the various calibration functions and their form. Specifically, we plotted the normalized linear terms (i.e., β1,i.i=1,,nη) of the Interaction energies and the bead size with respect to different values of the regularization term γ in Figs. 10(a) and 10(b), respectively. What we observe from these figures is that the parameters with larger values over the entire range of γ are the more important parameters, and functions, to achieve high calibration accuracy. Specifically, we find that ϵ2, ϵ3, and ϵ4 are highly influential to calibration accuracy. This is consistent with our observations in Fig. 5 where ϵ2, ϵ3 are approximated through a linear and RBF function, respectively. While ϵ4 has been approximated through a constant function, something that can be explained by noting that ϵ4 is highly correlated with ϵ3, as observed in Fig. 8, and its relatively small signal-to-noise ratio. Finally, we find that it is more important for the interaction energies to have more flexibility compared to the bead size. In summary, we find that we better calibration results can be achieved with the introduced semiparametric calibration method and also provides additional insight into the unknown calibration functions.

Fig. 10
Comparison of the normalized slopes for different regularization weights. The magnitude of the normalized slopes indicates the importance of different calibration functions with respect to calibration accuracy. In addition, the Interaction energies are most consequential to calibration accuracy compared to the bead size, an observation that is consistent with the quantified uncertainty obtained from the introduced semiparametric calibration approach. (a) Normalized slopes of the calibration functions for different regularization weights γ and Interaction energies (i.e., i ∈ 1,…, 7). (b) Normalized slopes of the calibration functions for different regularization weights γ and bead sizes (i.e., i ∈ 8,…, 14).
Fig. 10
Comparison of the normalized slopes for different regularization weights. The magnitude of the normalized slopes indicates the importance of different calibration functions with respect to calibration accuracy. In addition, the Interaction energies are most consequential to calibration accuracy compared to the bead size, an observation that is consistent with the quantified uncertainty obtained from the introduced semiparametric calibration approach. (a) Normalized slopes of the calibration functions for different regularization weights γ and Interaction energies (i.e., i ∈ 1,…, 7). (b) Normalized slopes of the calibration functions for different regularization weights γ and bead sizes (i.e., i ∈ 8,…, 14).
Close modal

While using Eq. (13) provides validation for the introduced semiparametric calibration approach, it also stipulates a path forward. Calibration results can be improved by providing a flexible class of calibration functions in Eq. (13) and by enabling autonomous selection. However, doing so is a nontrivial challenge, and as such, this will be of interest if future work.

5 Conclusions

In this paper, we presented a semiparametric functional calibration method for models with input-dependent parameters. Specifically, we introduced a modified Bayesian approximation scheme that uses the Hessian of a distance function to quantify the uncertainty in a high-dimensional candidate set of calibration functions. Subsequently, the quantified uncertainty was then used to define a more appropriate class of calibration functions with fewer parameters. The novelty of the presented method is the joint consideration of the following three modeling characteristics: (i) it allows the modeler to benefit from various sources of domain knowledge (e.g., the physics embedded in their simulations models, data from physical experiments, and expert knowledge), (ii) the modeler can greatly reduce the complexity of the calibration functions while maintaining strong predictive accuracy for high-dimensional problems (i.e., more than 10 calibration functions), and (iii) it can be considered as a semiparametric approach as it does not require the modeler to define the form of the calibration functions a priori. These characteristics make the introduced calibration method a powerful tool for modelers to explore physical phenomena at lengths and time scales unattainable with their high-fidelity data sources. This could have broad reaching implications on materials design, environmental modeling, and operations engineering. Finally, we have demonstrated the validity of the proposed method through an engineering problem where the number of calibration parameters was reduced by over 50% while only increasing the predictive error by 0.2%.

Future research in nonparametric functional calibration can be directed toward the consideration of a systematic bias between the low-fidelity and the high-fidelity data sources. Simulations often do not account for all the physics associated with the physical processes that they aim to model. Consequently, a systematic bias exists between the low-fidelity simulations and the high-fidelity physical experiments, necessitating the consideration of a bias function. However, the bias function introduces significant modeling flexibility that is likely to lead to identifiability issues. Consequently, a second direction for future research includes investigation into the calibration function identifiability. Available methods that address the identifiability issue involve stringent assumptions on the form of the bias function, necessitate an informative prior distribution, or require additional response variables that might not always be available. Third, the generalizability of the developed method could be tested by addressing challenges in other disciplines. For example, the modified approximate Bayesian computation scheme could be used to quantify the prediction uncertainty in the hyperparameters of deep learning models. Further research in Bayesian-based functional calibration will enable modelers to get a de-eper understanding of the systems that they are studying, and thus facilitate the design of products with improved utility.

Acknowledgment

Support from AFOSRFA9550-18-1-0381, U.S. Department of Commerce under Award No. 70 NANB19H005 and National Institute of Standards and Technology as part of the Center for Hierarchical Materials Design (CHiMaD). In addition, versions of Figs. 58 have originally been published in NPJ Computational Materials [57] under the creative commons attribution-noncommercial-sharealike license.

Funding Data

  • Air Force Office of Scientific Research (Award No. 70 NANB19H005; Funder ID: 10.13039/100000181).

  • National Institute of Standards and Technology (Award No. Chimad; Funder ID: 10.13039/100000161).

Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

Nomenclature

Acronyms
AA =

all atomistic

CG =

coarse-grained

DC =

degree of crosslinking

GP =

Gaussian process

MD =

molecular dynamics

RMSPE =

root mean squared percentage error

RBF =

radial basis function

UQ =

uncertainty quantification

Symbols
D =

distance metric function

d, p =

number of modeling variables and number of modeling parameters

E,Cov =

expectation and covariance function

f =

function and function value

k,K,w,c =

Vector of kernel function values, matrix of kernel function values, length scale parameter, and center values of RBF functions

m,M =

Eigenvector and column matrix of eigenvectors

n =

number of samples

N =

normal distribution

P =

probability

U,V,A =

a matrix of left singular vectors, a matrix of right singular vectors, and the matrix of singular values

X =

set of samples at the approximate Bayesian computation rejection boundary

γ =

smoothing parameter for regularization

ρ =

scaling parameter in the general calibration formulation

δ =

bias function

Σ =

covariance matrix

E =

space of admissible calibration functions

F =

space of possible function responses

η =

set of calibration functions

θ =

set of calibration parameters

ϑ =

reduced set of calibration parameters

ξ =

approximate Bayesian computation rejection value

μ,σ =

mean and standard deviation

Θ,Θ =

space of admissible calibration parameters and reduced dimensional calibration parameters

X,T =

space of admissible modeling variables and modeling parameters

ϵ,ε,e =

noise terms

Subscripts
l =

low-fidelity

h =

high-fidelity

L2 =

norm calculating a vectors Euclidean distance

Superscripts and Variable Accents
(k) =

index for the kth eigenvector

s,τ =

number of updates and total number of update of the sequential optimization scheme

=

simplified approximation

–1 =

inverse of a square matrix

* =

an optimized parameter or function

T =

vector transpose

̂  =

approximated parameters or function

˜  =

a variable arbitrary fixed at a specific value

Appendix: Discussion on the Introduced Method and Its Implementation

Concerning the proposed method's theoretical validity, as defined in Ref. [61], we continue our discussion by exploring some of its characteristics, limitations, and user scenarios in this section.

  1. In the introduced calibration method we approximate the first two statistical moments of the calibration functions. The uncertainty in these functions comes from the intrinsic data uncertainties as well as the interpolation uncertainty. By using the Hessian of the distance function at the global optimum we infer the covariance matrix of the calibration function parameters. However, in this process, we search for a rejection boundary that has a distance function value ξ=L(μ̂θ)/2 (Step 2c). This is a challenge also observed in approximate Bayesian computation [35] and only allows us to find a scaled version of the covariance matrix and as such also scales the range of the confidence bound of the approximated calibration functions. Consequently, the calibration method will only provide insight into the relative importance of calibration functions as opposed to absolute importance.

  2. An additional consideration is the potential multi-modality of the distance function. The introduced uncertainty quantification scheme infers the rejection boundary of the response around the global optimum of the calibration parameters and as such is unable to consider the local optimums isolated from the global optimum. Note that we used the term isolated, as the local optimum that is in the proximity of the global optimum will be explored by the calibration method. However, we believe that the influence of this limitation is marginal as the local minimum observed in our test problems has relatively small distance function values compared to the global optimum and as such would have a very small probability in a full Bayesian scheme (i.e., non approximate).

  3. The introduced method is most suitable for the scenario where we have a low-fidelity model and a high-fidelity model that given the true set of calibration functions has no systematic bias. When used for two simulation models this can be a reasonable assumption as they can be constructed to include the same set of physics. However, when we extend this to the calibration of a low-fidelity simulation model to a set of high-fidelity physical experiments, then the existence of a bias function needs to be considered. The reason for this is that those simulation models typically include only the most prevalent physics and as such are expected to have a systematic bias with physical experiments [12]. Modelers need to take care when using the introduced calibration method in such scenarios as the analysis will try to compensate for the missing physics by using the additional freedom provided by the candidate set of calibration functions.

  4. Another important consideration is that the calibration method requires the modeler to define a range for the calibration functions a priori. If this range is chosen to be too small then the true calibration function might extrapolate from this range and as such the identified calibration will provide only a little, or even misleading, information about the “true” form of the calibration functions. However, in most cases, this scenario is identified by finding a set of calibration functions that are at the calibration function bounds. In such a scenario, the modeler is advised to run additional low-fidelity simulations to increase the admissible range of the calibration functions.

  5. A final issue to discuss is identifiability [62,63], a phenomenon that occurs when there is either a non-unique set of optimal calibration functions or when there is a large range for the calibration function parameters that have similar calibration accuracy. The issue of identifiability can be observed in the introduced method when the variance (i.e., the diagonal elements of the calibration parameter matrix) is significantly larger than that of the other parameters (e.g., η6,η7 compared to η13 in the coarse-grained epoxy model as studied in Sec. 4.2). In such an event, the associate calibration function can in most cases be replaced by a constant function or linear function (e.g., η6(·) and η7(·) in Fig. 5). However, this does not imply that the “true” calibration function has a low-dimensional parametric form, we can only conclude that its form has little influence on the overall calibration accuracy compared to other calibration functions. We have found in the coarse-grained epoxy example that the identifiability of the model improved when calibrating all eight responses simultaneously as opposed to individually. This is an observation that makes intuitive sense as including multiple responses reduces the likelihood that there exists a non-unique solution [64].

The purpose of discussing the above issues is to bring attention to the various assumptions and considerations that underpin the presented calibration method. Moreover, knowledge of these issues informs how results should be interpreted and guides further scientific effort.

References

1.
Higdon
,
D.
,
Gattiker
,
J.
,
Williams
,
B.
, and
Rightley
,
M.
,
2008
, “
Computer Model Calibration Using High-Dimensional Output
,”
J. Am. Stat. Assoc.
,
103
(
482
), pp.
570
583
.10.1198/016214507000000888
2.
Leoni
,
N.
, and
Amon
,
C. H.
,
2000
, “
Bayesian Surrogates for Integrating Numerical, Analytical, and Experimental Data: Application to Inverse Heat Transfer in Wearable Computers
,”
IEEE Trans. Compon. Packaging Technol.
,
23
(
1
), pp.
23
32
.10.1109/6144.833038
3.
Higdon
,
D.
,
Kennedy
,
M.
,
Cavendish
,
J. C.
,
Cafeo
,
J. A.
, and
Ryne
,
R. D.
,
2004
, “
Combining Field Data and Computer Simulations for Calibration and Prediction
,”
SIAM J. Sci. Comput.
,
26
(
2
), pp.
448
466
.10.1137/S1064827503426693
4.
Gattiker
,
J.
,
Higdon
,
D.
,
Keller-McNulty
,
S.
,
McKay
,
M.
,
Moore
,
L.
, and
Williams
,
B.
,
2006
, “
Combining Experimental Data and Computer Simulations, With an Application to Flyerplate Experiments
,”
Bayesian Anal.
,
1
(
4
), pp.
765
792
.10.1214/06-BA125
5.
Xiong
,
Y.
,
Chen
,
W.
,
Tsui
,
K.-L.
, and
Apley
,
D. W.
,
2009
, “
A Better Understanding of Model Updating Strategies in Validating Engineering Models
,”
Comput. Methods Appl. Mech. Eng.
,
198
(
15–16
), pp.
1327
1337
.10.1016/j.cma.2008.11.023
6.
Li
,
Y. G.
,
2010
, “
Gas Turbine Performance and Health Status Estimation Using Adaptive Gas Path Analysis
,”
ASME J. Eng. Gas Turbines Power
,
132
(
4
), p.
041701
.10.1115/1.3159378
7.
Liu
,
X.
,
Mao
,
K.
,
Wang
,
X.
,
Wang
,
X.
, and
Wang
,
Y.
,
2020
, “
A Modified Quality Loss Model of Service Life Prediction for Products Via Wear Regularity
,”
Reliab. Eng. Syst. Safety
,
204
, p.
107187
.10.1016/j.ress.2020.107187
8.
Li
,
X.
,
Zhang
,
M.
,
Wang
,
Y.
,
Zhang
,
M.
,
Prasad
,
A.
,
Chen
,
W.
,
Schadler
,
L.
, and
Catherine Brinson
,
L.
,
2019
, “
Rethinking Interphase Representations for Modeling Viscoelastic Properties for Polymer Nanocomposites
,”
Materials
,
6
, p.
100277
.10.1016/j.mtla.2019.100277
9.
Xia
,
W.
,
Song
,
J.
,
Hansoge
,
N. K.
,
Phelan
,
F. R.
,
Keten
,
S.
, and
Douglas
,
J. F.
,
2018
, “
Energy Renormalization for Coarse-Graining the Dynamics of a Model Glass-Forming Liquid
,”
J. Phys. Chem. B
,
122
(
6
), pp.
2040
2045
.10.1021/acs.jpcb.8b00321
10.
Sacks
,
J.
,
Welch
,
W. J.
,
Mitchell
,
T. J.
, and
Wynn
,
H. P.
,
1989
, “
Design and Analysis of Computer Experiments
,”
Stat. Sci.
,
4
(
4
), pp.
409
423
.10.1214/ss/1177012413
11.
Han
,
G.
,
Santner
,
T. J.
, and
Rawlinson
,
J. J.
,
2009
, “
Simultaneous Determination of Tuning and Calibration Parameters for Computer Experiments
,”
Technometrics
,
51
(
4
), pp.
464
474
.10.1198/TECH.2009.08126
12.
Kennedy
,
M. C.
, and
O'Hagan
,
A.
,
2001
, “
Bayesian Calibration of Computer Models
,”
J. R. Stat. Soc.: Ser. B (Stat. Methodology)
,
63
(
3
), pp.
425
464
.10.1111/1467-9868.00294
13.
Qian
,
P. Z. G.
, and
Jeff Wu
,
C. F.
,
2008
, “
Bayesian Hierarchical Modeling for Integrating Low-Accuracy and High-Accuracy Experiments
,”
Technometrics
,
50
(
2
), pp.
192
204
.10.1198/004017008000000082
14.
Plumlee
,
M.
,
2017
, “
Bayesian Calibration of Inexact Computer Models
,”
J. Am. Stat. Assoc.
,
112
(
519
), pp.
1274
1285
.10.1080/01621459.2016.1211016
15.
Fugate
,
M.
,
Williams
,
B.
,
Higdon
,
D.
,
Hanson
,
K. M.
,
Gattiker
,
J.
,
Chen
,
S. R.
,
Unal
,
C.
,
Al
,
P.
,
Du
,
B.
, and
U-Nb
,
T.
,
2010
, “
Hierarchical Bayesian Analysis and the Preston-Tonks-Wallace Model
,” accessed Aug. 9, 2021, https://kmh-lanl.hansonhub.com/publications/larept05.pdf
16.
Ezzat
,
A. A.
,
Pourhabib
,
A.
, and
Ding
,
Y.
,
2018
, “
Sequential Design for Functional Calibration of Computer Models
,”
Technometrics
,
60
(
3
), pp.
286
296
.10.1080/00401706.2017.1377638
17.
Tao
,
S.
,
Apley
,
D. W.
,
Chen
,
W.
,
Garbo
,
A.
,
Pate
,
D. J.
, and
German
,
B. J.
,
2019
, “
Input Mapping for Model Calibration With Application to Wing Aerodynamics
,”
AIAA J.
,
57
(
7
), pp.
2734
2745
.10.2514/1.J057711
18.
Andrew Brown
,
D.
, and
Atamturktur
,
S.
,
2018
, “
Nonparametric Functional Calibration of Computer Models
,”
Statistica Sin.
,
28
(
2
), pp.
721
742
.10.5705/ss.202015.0344
19.
Plumlee
,
M.
,
Roshan Joseph
,
V.
, and
Yang
,
H.
,
2016
, “
Calibrating Functional Parameters in the Ion Channel Models of Cardiac Cells
,”
J. Am. Stat. Assoc.
,
111
(
514
), pp.
500
509
.10.1080/01621459.2015.1119695
20.
Pourhabib
,
A.
, and
Balasundaram
,
B.
,
2015
, “
Non-Isometric Curve to Surface Matching With Incomplete Data for Functional Calibration
,” Machine Learning, e-print arXiv:1508.01240.
21.
Farmanesh
,
B.
,
Pourhabib
,
A.
,
Balasundaram
,
B.
, and
Buchanan
,
A.
,
2021
, “
A Bayesian Framework for Functional Calibration of Expensive Computational Models Through Non-Isometric Matching
,”
IISE Trans.
,
53
(
3
), pp.
352
364
.10.1080/24725854.2020.1774688
22.
Jiang
,
Z.
,
Chen
,
S.
,
Apley
,
D. W.
, and
Chen
,
W.
,
2016
, “
Reduction of Epistemic Model Uncertainty in Simulation-Based Multidisciplinary Design
,”
ASME J. Mech. Des.
,
138
(
8
), p.
081403
.10.1115/1.4033918
23.
Lee
,
G.
,
Kim
,
W.
,
Oh
,
H.
,
Youn
,
B. D.
, and
Kim
,
N. H.
, oct
2019
, “
Review of Statistical Model Calibration and Validation–From the Perspective of Uncertainty Structures
,”
Struct. Multidiscip. Optim.
,
60
(
4
), pp.
1619
1644
.10.1007/s00158-019-02270-2
24.
Duan
,
K.
,
He
,
Y.
,
Li
,
Y.
,
Liu
,
J.
,
Zhang
,
J.
,
Hu
,
Y.
,
Lin
,
R.
,
Wang
,
X.
,
Deng
,
W.
, and
Li
,
L.
,
2019
, “
Machine-Learning Assisted Coarse-Grained Model for Epoxies Over Wide Ranges of Temperatures and Cross-Linking Degrees
,”
Mater. Des.
,
183
, p.
108130
.10.1016/j.matdes.2019.108130
25.
Taka
,
E.
,
Stein
,
S.
, and
Williamson
,
J. H.
,
2020
, “
Increasing Interpretability of Bayesian Probabilistic Programming Models Through Interactive Representations
,”
Front. Comput. Sci.
,
2
, p.
52
.10.3389/fcomp.2020.567344
26.
Picheny
,
V.
,
Ginsbourger
,
D.
,
Richet
,
Y.
, and
Caplin
,
G.
,
2013
, “
Quantile-Based Optimization of Noisy Computer Experiments With Tunable Precision
,”
Technometrics
,
55
(
1
), pp.
2
13
.10.1080/00401706.2012.707580
27.
Le Gratiet
,
L.
, and
Garnier
,
J.
,
2014
, “
Recursive co-Kriging Model for Design of Computer Experiments With Multiple Levels of Fidelity
,”
Int. J. Uncertainty Quantif.
,
4
(
5
), pp.
365
386
.10.1615/Int.J.UncertaintyQuantification.2014006914
28.
Kennedy
,
M.
, and
O'Hagan
,
A.
,
2000
, “
Predicting the Output From a Complex Computer Code When Fast Approximations Are Available
,”
Biometrika
,
87
(
1
), pp.
1
13
.10.1093/biomet/87.1.1
29.
Goh
,
J.
,
Bingham
,
D.
,
Holloway
,
J. P.
,
Grosskopf
,
M. J.
,
Kuranz
,
C. C.
, and
Rutter
,
E.
,
2013
, “
Prediction and Computer Model Calibration Using Outputs From Multifidelity Simulators
,”
Technometrics
,
55
(
4
), pp.
501
512
.10.1080/00401706.2013.838910
30.
Tuo
,
R.
, and
Jeff Wu
,
C. F.
,
2015
, “
Efficient Calibration for Imperfect Computer Models
,”
Ann. Stat.
,
43
(
6
), pp.
2331
2352
.10.1214/15-AOS1314
31.
Atamturktur
,
S.
,
Hegenderfer
,
J.
,
Williams
,
B.
,
Egeberg
,
M.
,
Lebensohn
,
R. A.
, and
Unal
,
C.
,
2015
, “
A Resource Allocation Framework for Experiment-Based Validation of Numerical Models
,”
Mech. Adv. Mater. Struct.
,
22
(
8
), pp.
641
654
.10.1080/15376494.2013.828819
32.
Schaback
,
R.
,
1999
, “
Native Hilbert Spaces for Radial Basis Functions i
,”
New Developments in Approximation Theory
,
M. W.
Müller
,
M. D.
Buhmann
,
D. H.
Mache
, and
M.
Felten
, eds.,
Birkhäuser Basel
,
Basel
, pp.
255
282
.
33.
Schölkopf
,
B.
,
Herbrich
,
R.
, and
Smola
,
A. J.
,
2001
, “
A Generalized Representer Theorem
,”
Computational Learning Theory
,
H.
David
and
W.
Bob
, eds.,
Springer Berlin Heidelberg
,
Berlin, Heidelberg
, pp.
416
426
.
34.
Rasmussen
,
C. E.
, and
Williams
,
C. K. I.
,
2005
,
Gaussian Processes for Machine Learning (Adaptive Computation and Machine Learning)
,
The MIT Press
, Cambridge, MA.
35.
Grazian
,
C.
, and
Fan
,
Y.
,
2020
, “
A Review of Approximate Bayesian Computation Methods Via Density Estimation: Inference for Simulator-Models
,”
WIREs Comput. Stat.
,
12
(
4
), p.
e1486
.10.1002/wics.1486
36.
Pavel
,
I.
,
Dmitrii
,
P.
,
Timur
,
G.
,
Dmitry
,
V.
, and
Wilson
,
A. G.
,
2018
, “
Averaging Weights Leads to Wider Optima and Better Generalization
,”
34th Conference on Uncertainty in Artificial Intelligence 2018,
Monterey, CA, Aug. 6, pp.
876
885
.
37.
Maddox
,
W. J.
,
Izmailov
,
P.
,
Garipov
,
T.
,
Vetrov
,
D. P.
, and
Wilson
,
A. G.
,
2019
, “
A Simple Baseline for Bayesian Uncertainty in Deep Learning
,” Proceedings of the 33rd International Conference on Neural Information Processing Systems,
Vancouver, BC, Canada, Dec. 8, pp.
13153
13164
.
38.
Plumlee
,
M.
, March
2014
, “
Fast Prediction of Deterministic Functions Using Sparse Grid Experimental Designs
,”
J. Am. Stat. Assoc.
,
109
(
508
), pp.
1581
1591
.10.1080/01621459.2014.900250
39.
Forrester
,
A. I.
, and
Keane
,
A. J.
,
2009
, “
Recent Advances in Surrogate-Based Optimization
,”
Prog. Aerosp. Sci.
,
45
(
1–3
), pp.
50
79
.10.1016/j.paerosci.2008.11.001
40.
Pritchard
,
J. K.
,
Seielstad
,
M. T.
,
Perez-Lezaun
,
A.
, and
Feldman
,
M. W.
,
1999
, “
Population Growth of Human Y Chromosomes: A Study of Y Chromosome Microsatellites
,”
Mol. Biol. Evol.
,
16
(
12
), pp.
1791
1798
.10.1093/oxfordjournals.molbev.a026091
41.
Marjoram
,
P.
,
Molitor
,
J.
,
Plagnol
,
V.
, and
Tavaré
,
S.
,
2003
, “
Markov Chain Monte Carlo Without Likelihoods
,”
Proc. Natl. Acad. Sci.
,
100
(
26
), pp.
15324
15328
.10.1073/pnas.0306899100
42.
Sisson
,
S. A.
,
Fan
,
Y.
, and
Beaumont
,
M. A.
,
2018
, “
Overview of Approximate Bayesian Computation
,”
arXiv:1802.09720
.10.48550/arXiv.1802.09720
43.
Arora
,
J.
,
2016
,
Introduction to Optimum Design
, 4nd ed.,
Elsevier
, London, UK.
44.
Bertsimas
,
D.
, and
Thiele
,
A.
,
2014
, “
Robust and Data-Driven Optimization: Modern Decision Making Under Uncertainty
,”
INFORMS Tutorials in Operations Research
, pp.
95
122
.
45.
Bertsimas
,
D.
,
Nohadani
,
O.
, and
Teo
,
K. M.
,
2010
, “
Robust Optimization for Unconstrained Simulation-Based Problems
,”
Oper. Res.
,
58
(
1
), pp.
161
178
.10.1287/opre.1090.0715
46.
Saltelli
,
A.
,
Annoni
,
P.
,
Azzini
,
I.
,
Campolongo
,
F.
,
Ratto
,
M.
, and
Tarantola
,
S.
,
2010
, “
Variance Based Sensitivity Analysis of Model Output. design and Estimator for the Total Sensitivity Index
,”
Comput. Physics Commun.
,
181
(
2
), pp.
259
270
.10.1016/j.cpc.2009.09.018
47.
Despotovic
,
M.
,
Nedic
,
V.
,
Despotovic
,
D.
, and
Despotovic
,
S.
,
2016
, “
Evaluation of Empirical Models for Predicting Monthly Mean Horizontal Diffuse Solar Radiation
,”
Renewable Sustainable Energy Rev.
,
56
, pp.
246
260
.10.1016/j.rser.2015.11.058
48.
Göçken
,
M.
,
Özçalici
,
M.
,
Boru
,
A.
, and
Dosdoğru
,
A. T.
,
2016
, “
Integrating Metaheuristics and Artificial Neural Networks for Improved Stock Price Prediction
,”
Expert Syst. Appl.
,
44
, pp.
320
331
.10.1016/j.eswa.2015.09.029
49.
Li
,
M.-F.
,
Tang
,
X.-P.
,
Wu
,
W.
, and
Liu
,
H.-B.
,
2013
, “
General Models for Estimating Daily Global Solar Radiation for Different Solar Radiation Zones in Mainland China
,”
Energy Convers. Manage.
,
70
, pp.
139
148
.10.1016/j.enconman.2013.03.004
50.
Schwarz
,
G.
,
1978
, “
Estimating the Dimension of a Model
,”
Ann. Stat.
,
6
(
2
), pp.
461
464
.10.1214/aos/1176344136
51.
Komarov
,
P. V.
,
Yu-Tsung
,
C.
,
Shih-Ming
,
C.
,
Khalatur
,
P. G.
, and
Reineker
,
P.
,
2007
, “
Highly Cross-Linked Epoxy Resins: An Atomistic Molecular Dynamics Simulation Combined With a Mapping/Reverse Mapping Procedure
,”
Macromolecules
,
40
(
22
), pp.
8104
8113
.10.1021/ma070702+
52.
Bandyopadhyay
,
A.
,
Valavala
,
P. K.
,
Clancy
,
T. C.
,
Wise
,
K. E.
, and
Odegard
,
G. M.
,
2011
, “
Molecular Modeling of Crosslinked Epoxy Polymers: The Effect of Crosslink Density on Thermomechanical Properties
,”
Polymer
,
52
(
11
), pp.
2445
2452
.10.1016/j.polymer.2011.03.052
53.
King
,
J. A.
,
Klimek
,
D. R.
,
Miskioglu
,
I.
, and
Odegard
,
G. M.
,
2013
, “
Mechanical Properties of Graphene Nanoplatelet/Epoxy Composites
,”
J. Appl. Polym. Sci.
,
128
(
6
), pp.
4217
4223
.10.1002/app.38645
54.
Wu
,
C.
, and
Xu
,
W.
,
2007
, “
Atomistic Simulation Study of Absorbed Water Influence on Structure and Properties of Crosslinked Epoxy Resin
,”
Polymer
,
48
(
18
), pp.
5440
5448
.10.1016/j.polymer.2007.06.038
55.
Meng
,
Z.
,
Bessa
,
M. A.
,
Xia
,
W.
,
Liu
,
W. K.
, and
Keten
,
S.
,
2016
, “
Predicting the Macroscopic Fracture Energy of Epoxy Resins From Atomistic Molecular Simulations
,”
Macromolecules
,
49
(
24
), pp.
9474
9483
.10.1021/acs.macromol.6b01508
56.
Li
,
C.
,
Browning
,
A. R.
,
Christensen
,
S.
, and
Strachan
,
A.
,
2012
, “
Atomistic Simulations on Multilayer Graphene Reinforced Epoxy Composites
,”
Compos. Part A: Appl. Sci. Manuf.
,
43
(
8
), pp.
1293
1300
.10.1016/j.compositesa.2012.02.015
57.
Giuntoli
,
A.
,
Hansoge
,
N. K.
,
van Beek
,
A.
,
Meng
,
Z.
,
Chen
,
W.
, and
Keten
,
S.
,
2021
, “
Systematic Coarse-Graining of Epoxy Resins With Machine Learning-Informed Energy Renormalization
,”
NPJ Comput. Mater.
,
7
(
1
), pp.
1
12
.10.1038/s41524-021-00634-1
58.
Sobol
,
I. M.
, May
1976
, “
Uniformly Distributed Sequences With an Additional Uniform Property
,”
USSR Comput. Math. Math. Phys.
,
16
(
5
), pp.
236
242
.10.1016/0041-5553(76)90154-3
59.
van Beek
,
A.
,
Tao
,
S.
, and
Chen
,
W.
,
2019
, “
Global Emulation Through Normative Decision Making and Thrifty Adaptive Batch Sampling
,”
ASME
Paper No. DETC2019-98223.10.1115/DETC2019-98223
60.
van Beek
,
A.
,
Tao
,
S.
,
Plumlee
,
M.
,
Apley
,
D. W.
, and
Chen
,
W.
,
2020
, “
Integration of Normative Decision-Making and Batch Sampling for Global Metamodeling
,”
ASME J. Mech. Des.
,
142
(
3
), p.
031114
.10.1115/1.4045601
61.
Pederson
,
K.
,
Emblemsvåg
,
J.
,
Bailey
,
R.
,
Allen
,
J. K.
, and
Mistree
,
F.
,
2000
, “
Validating Design Methods & Research: The Validation Square
,”
Design Engineering Technical Conferences
, Baltimore, MD, Sept. 10–14, pp.
1
13
.https://cecas.clemson.edu/cedar/wpcontent/uploads/2016/07/2-Mistree2000.pdf
62.
Arendt
,
P. D.
,
Apley
,
D. W.
, and
Chen
,
W.
,
2012
, “
Quantification of Model Uncertainty: Calibration, Model Discrepancy, and Identifiability
,”
ASME J. Mech. Des.
,
134
(
10
), p.
100908
.10.1115/1.4007390
63.
Jiang
,
Z.
,
Apley
,
D. W.
, and
Chen
,
W.
, January
2015
, “
Surrogate Preposterior Analyses for Predicting and Enhancing Identifiability in Model Calibration
,”
Int. J. Uncertainty Quantif.
,
5
(
4
), pp.
341
359
.10.1615/Int.J.UncertaintyQuantification.2015012627
64.
Arendt
,
P. D.
,
Apley
,
D. W.
,
Chen
,
W.
,
Lamb
,
D.
, and
Gorsich
,
D.
,
2012
, “
Improving Identifiability in Model Calibration Using Multiple Responses
,”
ASME J. Mech. Des.
,
134
(
10
), p.
100909
.10.1115/1.4007573