1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
|
#! /usr/bin/env python
import openturns as ot
import openturns.testing as ott
ot.TESTPREAMBLE()
ot.PlatformInfo.SetNumericalPrecision(3)
size = 100
# no observations
x = ot.Sample(size, 0)
g = ot.SymbolicFunction(
["a", "b", "c"],
[
"a + -1.0 * b + 1.0 * c",
"a + -0.6 * b + 0.36 * c",
"a + -0.2 * b + 0.04 * c",
"a + 0.2 * b + 0.04 * c",
"a + 0.6 * b + 0.36 * c",
"a + 1.0 * b + 1.0 * c",
],
)
outputDimension = g.getOutputDimension()
trueParameter = [2.8, 1.2, 0.5]
params = [0, 1, 2]
Theta1 = ot.Dirac(trueParameter[0])
Theta2 = ot.Dirac(trueParameter[1])
Theta3 = ot.Dirac(trueParameter[2])
inputRandomVector = ot.JointDistribution([Theta1, Theta2, Theta3])
model = ot.ParametricFunction(g, params, trueParameter)
inputSample = inputRandomVector.getSample(size)
y = g(inputSample)
outputObservationNoiseSigma = 0.05
meanNoise = ot.Point(outputDimension)
covarianceNoise = ot.Point(outputDimension, outputObservationNoiseSigma)
R = ot.IdentityMatrix(outputDimension)
observationOutputNoise = ot.Normal(meanNoise, covarianceNoise, R)
# Add noise
sampleNoise = observationOutputNoise.getSample(size)
y += sampleNoise
candidate = [1.0] * 3
bootstrapSizes = [0, 100]
for bootstrapSize in bootstrapSizes:
algo = ot.NonLinearLeastSquaresCalibration(model, x, y, candidate)
algo.setBootstrapSize(bootstrapSize)
algo.run()
# To avoid discrepance between the platforms with or without CMinpack
# Check MAP
calibrationResult = algo.getResult()
parameterMAP = calibrationResult.getParameterMAP()
print("(Auto) MAP=", repr(parameterMAP))
rtol = 1.0e-2
atol = 0.0
ott.assert_almost_equal(parameterMAP, trueParameter, rtol, atol)
multiStartSize = 10
algo.setOptimizationAlgorithm(
ot.MultiStart(
ot.TNC(),
ot.LowDiscrepancyExperiment(
ot.SobolSequence(),
ot.Normal(
candidate, ot.CovarianceMatrix(ot.Point(candidate).getDimension())
),
multiStartSize,
).generate(),
)
)
algo.run()
# To avoid discrepance between the platforms with or without CMinpack
calibrationResult = algo.getResult()
parameterMAP = calibrationResult.getParameterMAP()
print("(Multistart/TNC) MAP=", repr(parameterMAP))
rtol = 1.0e-2
atol = 0.0
ott.assert_almost_equal(parameterMAP, trueParameter, rtol, atol)
|
