1%feature("docstring") OT::LinearLeastSquaresCalibration
2"Linear least squares calibration algorithm.
3
4Available constructors:
5    LinearLeastSquaresCalibration(*model, inputObservations, outputObservations, candidate, methodName*)
6
7    LinearLeastSquaresCalibration(*modelObservations, gradientObservations, outputObservations, candidate, methodName*)
8
9Parameters
10----------
11model : :class:`~openturns.Function`
12    The parametric function to be calibrated.
13inputObservations : 2-d sequence of float
14    The sample of input observations.
15    Can have dimension 0 to specify no observations.
16outputObservations : 2-d sequence of float
17    The sample of output observations.
18candidate : sequence of float
19    The reference value of the parameter.
20methodName : str
21    The name of the least-squares method to use for the calibration.
22    By default, equal to *QR*. Possible values are *SVD*, *QR*, *Cholesky*.
23modelObservations : 2-d sequence of float
24    The sample of output values of the model.
25gradientObservations : 2-d sequence of float
26    The Jacobian matrix of the model with respect to the parameter.
27
28Notes
29-----
30LinearLeastSquaresCalibration is the minimum variance estimator of the parameter of a given model under the assumption that this parameter acts linearly in the model.
31
32The prior distribution of the parameter is a noninformative prior
33emulated using a flat :class:`~openturns.Normal` centered on the candidate and with a variance equal to SpecFunc.MaxScalar.
34
35The posterior distribution of the parameter is :class:`~openturns.Normal` and reflects the
36variability of the optimum parameter depending on the observation sample.
37The associated covariance matrix may be regularized depending on the value of the
38key `LinearLeastSquaresCalibration-Regularization` in the :class:`~openturns.ResourceMap`.
39Let us denote by :math:`s_1` the smallest singular value of the covariance matrix.
40The default value of the `LinearLeastSquaresCalibration-Regularization`, zero,
41ensures that the singular values of the covariance matrix are left unmodified.
42If this parameter is set to a nonzero, relatively small, value denoted by :math:`\epsilon`,
43then all singular values of the covariance matrix are increased by :math:`\epsilon s_1`.
44
45The resulting distribution of the output error is :class:`~openturns.Normal` with a zero mean
46and a diagonal covariance matrix computed from the residuals.
47The residuals are computed based on the linearization of the model,
48where the Jacobian matrix is evaluated at the candidate point.
49The diagonal of the covariance matrix of the output error
50is constant and is estimated with the unbiased variance estimator.
51
52See also
53--------
54GaussianLinearCalibration, NonLinearLeastSquaresCalibration, GaussianNonLinearCalibration
55
56Examples
57--------
58Calibrate a nonlinear model using linear least-squares:
59
60>>> import openturns as ot
61>>> ot.RandomGenerator.SetSeed(0)
62>>> m = 10
63>>> x = [[0.5 + i] for i in range(m)]
64>>> inVars = ['a', 'b', 'c', 'x']
65>>> formulas = ['a + b * exp(c * x)']
66>>> model = ot.SymbolicFunction(inVars, formulas)
67>>> p_ref = [2.8, 1.2, 0.5]
68>>> params = [0, 1, 2]
69>>> modelX = ot.ParametricFunction(model, params, p_ref)
70>>> y = modelX(x)
71>>> y += ot.Normal(0.0, 0.05).getSample(m)
72>>> candidate = [1.0]*3
73>>> method = 'SVD'
74>>> algo = ot.LinearLeastSquaresCalibration(modelX, x, y, candidate, method)
75>>> algo.run()
76>>> print(algo.getResult().getParameterMAP())
77[8.24019,0.0768046,0.992957]"
78
79// ---------------------------------------------------------------------
80
81%feature("docstring") OT::LinearLeastSquaresCalibration::getModelObservations
82"Accessor to the model evaluation at the candidate.
83
84Returns
85-------
86modelObservation : :class:`~openturns.Sample`
87    Evaluation of the model at the candidate point."
88
89// ---------------------------------------------------------------------
90
91%feature("docstring") OT::LinearLeastSquaresCalibration::getGradientObservations
92"Accessor to the model gradient at the candidate.
93
94Returns
95-------
96gradientObservation : :class:`~openturns.Matrix`
97    Gradient of the model at the candidate point."
98
99// ---------------------------------------------------------------------
100
101%feature("docstring") OT::LinearLeastSquaresCalibration::getCandidate
102"Accessor to the parameter candidate.
103
104Returns
105-------
106candidate : :class:`~openturns.Point`
107    Parameter candidate."
108
109// ---------------------------------------------------------------------
110
111%feature("docstring") OT::LinearLeastSquaresCalibration::getMethodName
112"Accessor to the name of least-squares method used for the resolution.
113
114Returns
115-------
116name : :class:`~openturns.String`
117    Name of least-squares method used for the resolution."
118
119