An overview of the Scidoe toolbox.
The goal of this toolbox is to provide design of experiments (DOE) techniques, along with functions for model building. The provided functions are designed to be compatible with Matlab.
The toolbox is based on macros, which are vectorized for performances.
In this tutorial example, we combine the DOE features of
the Scidoe toolbox with the regression features in the Stixbox
toolbox.
More precisely, we consider a back-box function
myfunction
with 2 inputs and 1 output and try to
create a response surface which approximates this function.
We consider here a Central Composite Design (CCD) of experiments,
computed with the scidoe_ccdesign
function.
Then we evaluate the function at the design points and
use the regres
function in order to
fit a nonlinear model.
We begin by creating the design of experiments with the
scidoe_ccdesign
function.
This function produces points in the [-1,1] multidimensional
interval.
A CCD is composed of a factorial block design,
a block of axial (star) points and a matrix of center points (zeros).
In general, this design is appropriate for physical experiments
where a measure at a particular input point does not produce
the same output, because there are measure errors.
Hence, this design is associated with several points at the center
of the domain.
The previous script produces the following output.
-->x = scidoe_ccdesign(2) x = - 1. - 1. - 1. 1. 1. - 1. 1. 1. 0. 0. 0. 0. 0. 0. 0. 0. - 1.4142136 0. 1.4142136 0. 0. - 1.4142136 0. 1.4142136 0. 0. 0. 0. 0. 0. 0. 0.
The number of experiments corresponds to the number of rows : there are 16 experiments in this design.
-->nbexp=size(x,"r") nbexp = 16.
In order to see what points have been selected, we create a plot with the following statements.
The previous script produces the following output.
In this design, there are 8 points at the center of the domain and we might as well get only one point. This can be done with the "center" option, which specifies the number of central points. The following statement produces only one central point in the factorial design (the "1") and zero central point in the axial design (the "0").
With these settings, the number of experiments reduces to 9.
-->x = scidoe_ccdesign(2,"center",[1 0]) x = - 1. - 1. - 1. 1. 1. - 1. 1. 1. 0. 0. - 1.2649111 0. 1.2649111 0. 0. - 1.2649111 0. 1.2649111 -->nbexp=size(x,"r") nbexp = 9.
We then define our "black-box" function. This function involves a constant, two terms linear in x(1) and x(2) and a nonlinear term x(1)*x(2). The function is vectorized, so that we can evaluate the whole DOE in just one statement.
We add some amount of random noise by adding to the output random normal outcomes with mean zero and standard deviation 0.1.
We now have a column matrix of randomized function values.
-->y y = 8.0614563 - 9.8526572 4.1690109 5.9797968 1.8529469 - 1.8283189 5.9630483 7.0722875 - 3.0547834
The matrix y contains data which may have been either measured, or computed by the means of a complicated, external, software, such as a 3D finite element simulation, for example. In other words, although the current example is simple and analytic, it can be applied to much more complicated situations where the function to approximate is unknown.
Then we create the regression model, by using the
regres
function.
We are searching for the parameters B of the model:
y=B(1)+B(2)*x(1)+B(3)*x(2)+B(4)*x(1)*x(2)
This model is linear in the parameter B.
The regres
function solves a linear least squares
problems so that the squared differences between the
model and the data are minimized.
B = regres(y,[ones(y),x,x(:,1).*x(:,2)]) | ![]() | ![]() |
The previous script produces:
-->B = regres(y,[ones(y),x,x(:,1).*x(:,2)]) B = 2.0403097 3.0271382 - 4.0157183 4.9312249
We see that the estimated coefficients are close to the exact
ones.
In fact, if we do not introduce a random noise, then the
computation is been almost exact in the sense that the
coefficients are extremely close to the exact ones.
We now define the function mymodel
which uses the coefficients B to compute the
output of the regression model y2
.
function y=mymodel(x, B) y=B(1)+B(2)*x(:,1)+B(3)*x(:,2)+B(4)*x(:,1).*x(:,2) endfunction y2=mymodel(x,B); | ![]() | ![]() |
In order to see how the outputs from the model fit the outputs from the blackbox, we make the following graphics, where we plot the data versus the model outputs.
The previous script produces the following output.