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
87
88
89
90
91
92
93
94
95
|
"""
Build and validate a linear model
=================================
"""
# %%
# In this example, we build a linear regression model and validate it numerically and graphically.
#
# The linear model links a scalar variable :math:`Y` and to an n-dimensional one :math:`\underline{X} = (X_i)_{i \leq n}`, as follows:
#
# .. math::
# \tilde{Y} = a_0 + \sum_{i=1}^n a_i X_i + \varepsilon
#
# where :math:`\varepsilon` is the residual, supposed to follow :math:`\mathcal{N}(0.0, 1.0)`.
#
# The linear model may be validated graphically if :math:`\underline{X}` is of dimension 1, by drawing on the same graph the cloud :math:`(X_i, Y_i)`.
#
# The linear model can also be validated numerically with several tests:
#
# - `LinearModelFisher`: tests the nullity of the regression linear model coefficients (Fisher distribution used),
# - `LinearModelResidualMean`: tests, under the hypothesis of a Gaussian sample, if the mean of the residual is equal to zero.
# It is based on the Student test (equality of mean for two Gaussian samples).
#
#
# The hypothesis on the residuals (centered Gaussian distribution) may be validated:
#
# - graphically if :math:`\underline{X}` is of dimension 1, by drawing the residual couples (:math:`\varepsilon_i, \varepsilon_{i+1}`),
# where the residual :math:`\varepsilon_i` is evaluated on the samples :math:`(X, Y)`.
# - numerically with the `LinearModelResidualMean` test which tests, under the hypothesis of a Gaussian sample, if the mean of the residual is equal to zero.
# It is based on the Student test (equality of mean for two Gaussian samples).
#
# %%
import openturns as ot
import openturns.viewer as otv
# %%
# Generate `X, Y` samples
N = 1000
Xsample = ot.Triangular(1.0, 5.0, 10.0).getSample(N)
Ysample = Xsample * 3.0 + ot.Normal(0.5, 1.0).getSample(N)
# %%
# Generate a particular scalar sampleX
particularXSample = ot.Triangular(1.0, 5.0, 10.0).getSample(N)
# %%
# Create the linear model from `Y, X` samples
result = ot.LinearModelAlgorithm(Xsample, Ysample).getResult()
# Get the coefficients ai
print("coefficients of the linear regression model = ", result.getCoefficients())
# Get the confidence intervals of the `ai` coefficients
print(
"confidence intervals of the coefficients = ",
ot.LinearModelAnalysis(result).getCoefficientsConfidenceInterval(0.9),
)
# %%
# Validate the model with a visual test
graph = ot.VisualTest.DrawLinearModel(Xsample, Ysample, result)
view = otv.View(graph)
# %%
# Draw the graph of the residual values
graph = ot.VisualTest.DrawLinearModelResidual(Xsample, Ysample, result)
view = otv.View(graph)
# %%
# Check the nullity of the regression linear model coefficients
resultLinearModelFisher = ot.LinearModelTest.LinearModelFisher(
Xsample, Ysample, result, 0.10
)
print("Test Success ? ", resultLinearModelFisher.getBinaryQualityMeasure())
print("p-value of the LinearModelFisher Test = ", resultLinearModelFisher.getPValue())
print("p-value threshold = ", resultLinearModelFisher.getThreshold())
# %%
# Check, under the hypothesis of a Gaussian sample, if the mean of the residuals is equal to zero
resultLinearModelResidualMean = ot.LinearModelTest.LinearModelResidualMean(
Xsample, Ysample, result, 0.10
)
print("Test Success ? ", resultLinearModelResidualMean.getBinaryQualityMeasure())
print(
"p-value of the LinearModelResidualMean Test = ",
resultLinearModelResidualMean.getPValue(),
)
print("p-value threshold = ", resultLinearModelResidualMean.getThreshold())
# %%
otv.View.ShowAll()
|
