Post on 30-Aug-2019
Beispiel: Lineare Regression mit
Mathematica
Daten aus einem csv-File einlesen
In[1]:= data = Import@"�Users�ihn�Documents�Teaching�VP-Leitung�DataAnalysis�2013NewProgram�
Lectures given in 2013�6. Lecture�daten.csv"D
Out[1]= 881, 4.52318<, 82, 6.07091<, 83, 9.60907<, 84, 14.1273<,
85, 15.3875<, 86, 18.9468<, 87, 20.6456<, 88, 23.6468<<
Mit Datenlisten arbeiten
Länge der Datenliste (Zahl der Datenpunkte) ermitteln:
In[2]:= n = Length@dataD
Out[2]= 8
Auf einzelne Elemente (Datenpunkte) der Datenliste zugreifen:
In[3]:= data@@3DD
Out[3]= 83, 9.60907<
Auf Bereiche in der Datenliste zugreifen:
In[4]:= data@@Range@2, 4DDD
Out[4]= 882, 6.07091<, 83, 9.60907<, 84, 14.1273<<
Nur die x-Werte der Datenpunkte 2 bis 4:
In[5]:= data@@Range@2, 4D, 1DD
Out[5]= 82, 3, 4<
Nur die y-Werte der Datenpunkte 2 bis 4:
In[6]:= data@@Range@2, 4D, 2DD
Out[6]= 86.07091, 9.60907, 14.1273<
Datenliste graphisch darstellen
In[7]:= g1 = ListPlot@data, PlotRange ® 880, 9<, 80, 25<<,
PlotMarkers ® Automatic, LabelStyle ® 24, AxesLabel ® 8"x", "y"<D
Out[7]=
æ
æ
æ
æ
æ
æ
æ
æ
0 2 4 6 8x0
5
10
15
20
25y
Statistik der Datenliste
Mittelwert der x-Werte:
In[8]:= Meanx = Mean@data@@Range@1, nD, 1DDD
Out[8]=
9
2
Mittelwert der y-Werte:
In[9]:= Meany = Mean@data@@Range@1, nD, 2DDD
Out[9]= 14.1196
Varianz der x-Werte:
In[10]:= Varx = CentralMoment@data@@Range@1, nD, 1DD, 2D
Out[10]=
21
4
Varianz der y - Werte :
In[11]:= Vary = CentralMoment@data@@Range@1, nD, 2DD, 2D
Out[11]= 41.9354
Empirischer Korrelationskoeffizient:
In[12]:= Ρ =
CentralMoment@data, 81, 1<D
Varx Vary
Out[12]= 0.994133
2 LinearRegressionMathematicaForStudents.nb
Lineare Regression: Schätzwerte und ihre Unsicherheiten (gemäss Vorlesung)
In[13]:= A0 =
Vary
Varx
Ρ
Out[13]= 2.80967
In[14]:= B0 = Meany
Out[14]= 14.1196
In[15]:= C0 = Vary I1 - Ρ2M
Out[15]= 0.490655
In[16]:= ΣA =
C0
Varx Hn - 2L
Out[16]= 0.124805
In[17]:= ΣB =
C0
Hn - 2L
Out[17]= 0.285965
Fitkurve graphisch darstellen
In[18]:= g2 = Plot@A0 Hx - MeanxL + B0, 8x, 0, 8<,
PlotStyle ® Red, AxesLabel ® 8"x", "y"<, LabelStyle ® 24D
Out[18]=
2 4 6 8x
5
10
15
20
y
LinearRegressionMathematicaForStudents.nb 3
Daten mit Fitkurve graphisch darstellen
In[19]:= Show@g1, g2D
Out[19]=
æ
æ
æ
æ
æ
æ
æ
æ
0 2 4 6 8x0
5
10
15
20
25y
Variablen löschen
In[20]:= n =.;
Mit Wahrscheinlichkeitsdichteverteilungen arbeiten
Formel für bestimmte Verteilungen, die wir bereits kennen:
In[21]:= PDF@NormalDistribution@Μ, ΣD, xD
Out[21]=
ã-
Hx-ΜL2
2 Σ2
2 Π Σ
In[22]:= PDF@GammaDistribution@Α, ΒD, xD
Out[22]=
ã-
x
Β x-1+Α Β-Α
Gamma@ΑDx > 0
0 True
In[23]:= PDF@PoissonDistribution@ΜD, xD
Out[23]=
ã-Μ Μx
x!x ³ 0
0 True
In[24]:= PDF@BinomialDistribution@n, pD, kD
Out[24]=
H1 - pL-k+np
kBinomial@n, kD 0 £ k £ n
0 True
In[25]:= PDF@NegativeBinomialDistribution@n, pD, kD
Out[25]=
H1 - pLkp
nBinomial@-1 + k + n, -1 + nD k ³ 0
0 True
4 LinearRegressionMathematicaForStudents.nb
In[26]:= PDF@StudentTDistribution@Μ, Σ, nD, xD
Out[26]=
n
n+Hx-ΜL2
Σ2
1+n
2
n Σ BetaB n
2,
1
2F
In[27]:= PDF@InverseGammaDistribution@Α, ΒD, xD
Out[27]=
ã-
Β
x JΒ
xN
Α
x Gamma@ΑDx > 0
0 True
Mittelwerte von Verteilungen:
In[28]:= Mean@StudentTDistribution@Μ, Σ, nDD
Out[28]= ¶ Μ n > 1
Indeterminate True
Varianz von Verteilungen :
In[29]:= Variance@StudentTDistribution@Μ, Σ, nDD
Out[29]=
n Σ2
-2+nn > 2
Indeterminate True
Wahrscheinlichkeitsverteilungen graphisch darstellen:
In[30]:= Plot@PDF@StudentTDistribution@1, 2, 6D, xD,
8x, -13, 15<, AxesLabel ® 8"x", "pdfHxL"<, LabelStyle ® 24D
Out[30]=
-10 -5 5 10 15x
0.05
0.10
0.15
pdfHxL
Graphen animieren mit "Manipulate":
LinearRegressionMathematicaForStudents.nb 5
In[31]:= Manipulate@Plot@PDF@StudentTDistribution@Μ, Σ, nD, xD, 8x, -13, 15<,
AxesLabel ® 8"x", "pdfHxL"<, LabelStyle ® 24, PlotRange ® 88-15, 15<, 80, 1<<D,
88Μ, 0.5<, -10, 10<, 88Σ, 1<, 0.1, 2<, 88n, 1<, 1, 10, 1<D
Out[31]=
Μ
Σ
n
-15 -10 -5 0 5 10 15x
0.2
0.4
0.6
0.8
1.0pdfHxL
Die Verteilung für Σ2
In[32]:= n = Length@dataD
Out[32]= 8
In[33]:= PlotBPDFBInverseGammaDistributionBn
2
,
n C0
2
F, xF, 8x, 0, 4<,
PlotRange ® All, AxesLabel ® 9"Σ2", "pdfHΣ
2Èdata,M,IL"=,
LabelStyle ® 24, PlotStyle ® Thickness@0.01D, ImageSize ® 450F
Out[33]=
1 2 3 4Σ
2
0.5
1.0
1.5
pdfHΣ2Èdata,M,IL
6 LinearRegressionMathematicaForStudents.nb
Mittelwert für Σ2:
In[34]:=
n C0
n - 2
Out[34]= 0.654207
Standardabweichung für Σ2:
In[35]:=
n C0
n - 2
2
n - 4
Out[35]= 0.462594
Eine Abkürzung zum Linear Model Fitting in Mathematica
Mit dem folgenden “high-level” Befehl kann man eine lineare Regressionsgerade an Daten fitten:
In[36]:= model = LinearModelFit@data, x, xD
Out[36]= FittedModelB 1.47614 + 2.80967 x F
Fitkurve extrahieren:
In[37]:= model@"BestFit"D
Out[37]= 1.47614 + 2.80967 x
Fitkurve darstellen:
In[38]:= Show@ListPlot@data, PlotMarkers ® AutomaticD,
Plot@model@"BestFit"D, 8x, 0, 9<, PlotStyle ® RedD, PlotRange ® 880, 9<, 80, 25<<,
AxesOrigin ® 80, 0<, AxesLabel ® 8"x", "y"<, LabelStyle ® 24D
Out[38]=
æ
æ
æ
æ
æ
æ
æ
æ
2 4 6 8x
5
10
15
20
25
y
Die Fitparameter und ihre Fehler lassen sich so extrahieren :
In[39]:= model@"ParameterTable"D
Out[39]=
Estimate Standard Error t-Statistic P-Value
1 1.47614 0.630236 2.3422 0.0576706
x 2.80967 0.124805 22.5124 5.0275 ´ 10-7
Aber Vorsicht: der Fit mit der Funktion f(x) = a x + b führt zu korrelierten Parametern a und b. Die
Angabe des “Standard Error” reicht daher nicht aus, vielmehr muss eine Korrelationsmatrix
angegeben werden. Dieses Problem entsteht nicht beim Fitten mit f(x) = A (x - x) + B, da in diesem
Fall die Parameter A und B unkorreliert sind!
LinearRegressionMathematicaForStudents.nb 7
Aber Vorsicht: der Fit mit der Funktion f(x) = a x + b führt zu korrelierten Parametern a und b. Die
Angabe des “Standard Error” reicht daher nicht aus, vielmehr muss eine Korrelationsmatrix
angegeben werden. Dieses Problem entsteht nicht beim Fitten mit f(x) = A (x - x) + B, da in diesem
Fall die Parameter A und B unkorreliert sind!
Hier wäre ein Weg drumherum : wir verschieben die Daten “von Hand” um den Mittelwert von x.
In[40]:= datashifted = Table@8data@@j, 1DD - Meanx, data@@j, 2DD<, 8j, 1, n<D
Out[40]= ::-7
2
, 4.52318>, :-5
2
, 6.07091>, :-3
2
, 9.60907>, :-1
2
, 14.1273>,
:1
2
, 15.3875>, :3
2
, 18.9468>, :5
2
, 20.6456>, :7
2
, 23.6468>>
In[41]:= model = LinearModelFit@datashifted, x, xD
Out[41]= FittedModelB 14.1196 + 2.80967 x F
Jetzt stimmen die Fitparameter mit unseren Ergebnissen oben überein, und auch der “Standard
Error” gibt unsere Ergebnisse von oben genau wieder.
In[42]:= model@"ParameterTable"D
Out[42]=
Estimate Standard Error t-Statistic P-Value
1 14.1196 0.285965 49.3755 4.62844 ´ 10-9
x 2.80967 0.124805 22.5124 5.0275 ´ 10-7
Fazit: benutze keine Fitroutinen irgendeiner Software, von der Du nicht genau weisst, was
sie macht und die Du nicht an einfachsten Beispielen getestet hast! Die Zahl unkritischer
Anwender von Statistiksoftware ist bereits gross genug, die Zahl der grob falschen Ergeb-
nisse ist unermesslich.
8 LinearRegressionMathematicaForStudents.nb