(*^ ::[ frontEndVersion = "Microsoft Windows Mathematica Notebook Front End Version 2.2"; microsoftWindowsStandardFontEncoding; fontset = title, "Helv", 24, L0, center, nohscroll, bold; fontset = subtitle, "Helv", 18, L0, center, nohscroll, bold; fontset = subsubtitle, "Helv", 14, L0, center, nohscroll, bold; fontset = section, "Helv", 14, L0, bold, grayBox; fontset = subsection, "Helv", 12, L0, bold, blackBox; fontset = subsubsection, "Helv", 10, L0, bold, whiteBox; fontset = text, "Helv", 12, L0; fontset = smalltext, "Helv", 10, L0; fontset = input, "Courier New", 12, L0, nowordwrap, bold; fontset = output, "Courier New", 12, L0, nowordwrap; fontset = message, "Courier New", 10, L0, nowordwrap, R65280; fontset = print, "Courier New", 10, L0, nowordwrap; fontset = info, "Courier New", 10, L0, nowordwrap; fontset = postscript, "Courier New", 8, L0, nowordwrap; fontset = name, "Helv", 10, L0, nohscroll, italic, B65280; fontset = header, "Helv", 18, L0, nohscroll, bold; fontset = footer, "Helv", 18, L0, center, nohscroll, bold; fontset = help, "Helv", 10, L0, nohscroll; fontset = clipboard, "Helv", 12, L0, nohscroll; fontset = completions, "Helv", 12, L0, nowordwrap, nohscroll; fontset = graphics, "Courier New", 10, L0, nowordwrap, nohscroll; fontset = special1, "Arial", 12, L0, center, nowordwrap, nohscroll, bold; fontset = special2, "Helvetica", 12, L0, center, nowordwrap, nohscroll, bold; fontset = special3, "Helv", 12, L0, right, nowordwrap, nohscroll; fontset = special4, "Helv", 12, L0, nowordwrap, nohscroll; fontset = special5, "Helv", 12, L0, nowordwrap, nohscroll; fontset = leftheader, "Helv", 12, L0, nowordwrap, nohscroll; fontset = leftfooter, "Helv", 12, L0, nowordwrap, nohscroll; fontset = reserved1, "Courier New", 10, L0, nowordwrap, nohscroll;] :[font = title; inactive; nohscroll; center; backColorRed = 65280; backColorGreen = 65280; backColorBlue = 65280; fontColorRed = 32768; fontColorGreen = 0; fontColorBlue = 0; bold; fontName = ""; fontSize = 24; ] Statistik mit Mathematica :[font = subsection; inactive; startGroup; Cclosed; ] Copyright :[font = text; inactive; endGroup; ] Copyright 1994, Claudia Funke TU-Berlin, FB 13, Fachgebiet Ökonometrie und Statistik Dieses Notebook darf ausschließlich als Unterrichtsmaterial und für private Zwecke verwendet und nicht ohne Zustimmung der Autorin verändert werden. :[font = section; inactive; startGroup; Cclosed; ] Ausgangsituation :[font = text; inactive; ] Nach den jüngsten Veröffentlichungen über die Schädlichkeit des Rauchens macht sich die Firmenleitung einer größeren Gesellschaft Gedanken, wie sie die bei ihr beschäftigten Raucher und Raucherinnen positiv beeinflußen kann. Dazu läßt sie sich von ihrer statistischen Abteilung mittels einer repräsentativen Stichprobe das Rauchverhalten ihrer Beschäftigten ermitteln. :[font = text; inactive; ] Diese wurden dementsprechend in fünf von der Firmenleitung vorgegebenen Beschäftigungsgruppen aufgeteilt befragt. Die Werte wurden anschließend in zuvor definierte Raucherklassen eingeordnet: nicht = 0 Zigaretten / Tag leicht = 1 - 3 Zigaretten / Tag mittel = 4 - 19 Zigaretten / Tag stark = 19 - 35 Zigaretten / Tag :[font = text; inactive; ] Es ergab sich danach folgende Gesamtübersicht: :[font = text; inactive; output; endGroup; backColorRed = 65280; backColorGreen = 65280; backColorBlue = 65280; fontColorRed = 0; fontColorGreen = 0; fontColorBlue = 0; plain; fontName = "Arial"; fontSize = 12; ] RAUCHERTYPUS nicht leicht mittel stark --------------------------- Jüng.Arb. 36 48 66 26 Ält. Arb. 50 20 24 0 Jüng. Ang. 0 6 14 0 Ält. Ang. 8 4 6 4 Management 20 12 14 4 ------------------------------------------------------------------------------------ Die MitarbeiterInnen der statistischen Abteilung haben nun die Aufgabe, diese Daten statisch aufzubereiten. :[font = section; inactive; startGroup; ] Deskriptive Statistik :[font = subsubsection; inactive; startGroup; Cclosed; ] Einlesen und Darstellen eines Datensatzes :[font = text; inactive; ] Die als ASCII - File vorliegenden numerischen Daten werden mit :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] ReadList[ "file", Number, RecordLists-> True ] :[font = text; inactive; ] in Mathematica eingelesen. Diese spezielle Form der Funktion ReadList[ ] liest eine Folge von Zahlen (Number) ein, die durch ein Leerzeichen voneinander getrennt sind. Die Option RecordLists -> True bewirkt, daß ein als zweidimensionale Matrix in genau diesem Format eingelesen wird. :[font = input; nowordwrap; ] xKlass = ReadList["rauch.dat", Number, RecordLists -> True]; MatrixForm[xKlass] :[font = text; inactive; ] Im Gegensatz zur Matrixdarstellung läßt die Darstellung einer mehrdimensionalen Liste mit :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] TableForm[ Liste ] :[font = text; inactive; ] mehr Gestaltungsmöglichkeiten zu, z.B. Bezeichnungen für die Zeilen und Spalten einer Matrix: :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] TableForm[ Liste, TableHeadings -> {"Text1", "Text2", ... } ] :[font = input; endGroup; nowordwrap; ] TableForm[xKlass, TableHeadings -> \ { {"jung. Arb.","alt. Arb.","jung. Angst.", "alt. Angst.","Manag."}, {"nicht","leicht","mittel", "stark"} }, \ TableDirections -> {Column, Row}, \ TableAlignments -> {Center,Right} \ ] :[font = subsubsection; inactive; startGroup; Cclosed; ] Analyse des Datensatzes :[font = text; inactive; ] Einen repräsentativen Überblick über die Situation der Nichtraucher und Nichtraucherinnen erhält man, indem man die Daten für die Gesamt- belegschaft betrachtet. :[font = subsubsection; inactive; startGroup; ] Umstrukturierung des Datensatzes :[font = text; inactive; ] Ddie Daten der gesamten Belegschaft (für die einzelnen Rauchertypen) ergeben sich aus der Randhäufigkeit der Spalten: :[font = special2; inactive; output; nowordwrap; nohscroll; center; ] n h_j = Sum h_ij i=1 :[font = text; inactive; ] mit h_ij = Häufigkeit einer Merkmalkombination in der i - ten Zeile und der j - ten Spalte. :[font = text; inactive; ] Die Summen aller Spalten der Datenmatrix erhält man in Mathematica am elegantesten mit einer Funktion, die die benötigte Operation auf mehrere Listen anwendbar macht. Mit :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] Apply[ Plus, Liste ] :[font = text; inactive; ] oder :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] Plus @@ Liste :[font = text; inactive; ] werden die Elemente jeder Liste in list aufaddiert. :[font = input; nowordwrap; ] belegKlass = Apply[Plus, xKlass] :[font = text; inactive; ] Der Umfang der Stichprobe wiederum ist die Summe der Randhäufigkeiten. :[font = special2; inactive; output; nowordwrap; nohscroll; center; ] k N = Sum h_j j=1 :[font = input; nowordwrap; ] spumfang = Sum[belegKlass[[i]], {i, Length[belegKlass]} ] :[font = text; inactive; ] Die neue Zeile belegKlass läßt sich mit der Funktion AppendRows[ ] aus dem Package LinearAlgebra`MatrixManipulation` in die Datenmatrix X einfügen. :[font = input; nowordwrap; ] < \ { {"jung. Arb.","alt. Arb.","jung. Angst.", "alt. Angst.","Manag.","Gesamt"}, \ {"nicht","leicht","mittel","stark"} }, TableDirections -> {Column, Row}, \ TableAlignments -> {Center,Right} \ ] :[font = subsubsection; inactive; startGroup; ] Graphische Darstellungsmöglichkeiten :[font = text; inactive; ] Die Werte für die Gesamtbelegschaft lassen sich anschaulich in einem Balkendiagramm (absolute Häufigkeit) und in einem Tortendiagramm (relative Häufigkeit) darstellen. Die entsprechenden Funktionen sind im Package Graphics`Graphics` enthalten. :[font = input; nowordwrap; ] < {xticks, yticks} ermöglicht z.B. die Bennenung einzelner Balken, die Option PlotLabel -> "String" die Beschriftung der Graphik. :[font = input; nowordwrap; ] barpic = BarChart[belegKlass, Ticks -> \ { {{1,"nicht"}, {2,"leicht"}, {3,"mittel"}, {4,"stark"}}, \ Automatic}, PlotLabel -> "Rauchertypen" \ ] :[font = text; inactive; ] Mit der Funktion InputForm [ ] kann man hier sehr schön sehen, wie Mathematica solch eine Graphik aus einzelnen Graphikelementen aufbaut: :[font = input; nowordwrap; ] InputForm[barpic] :[font = input; nowordwrap; ] Mit :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] PieChart[ {y1, y2, ... } ] :[font = text; inactive; ] wird ein Tortendiagramm für die Raucherklassen bezüglich der Gesamt- belegschaft erzeugt. :[font = input; endGroup; endGroup; nowordwrap; ] tortpic = PieChart[{ {belegKlass[[1]], "nicht"}, {belegKlass[[2]], "leicht"}, {belegKlass[[3]], "mittel"}, {belegKlass[[4]], "stark"} } ] :[font = subsubsection; inactive; startGroup; Cclosed; ] Statistische Kenngrößen :[font = text; inactive; ] Der erste graphische Eindruck soll nun durch statistische Kenngrößen quantifiziert werden. :[font = text; inactive; ] Doch dazu werden statt der klassierten Daten die ursprünglichen Stichprobenwerte verwandt (Zigaretten / Tag). Dieser Datensatz wird wiederum mit ReadList[ ] eingelesen: :[font = input; nowordwrap; ] x = ReadList["rauchorg.dat", Number, RecordLists -> True] :[font = input; nowordwrap; ] Dimensions[x] :[font = text; inactive; ] Die Werte für die Gesamtbelegschaft lauten damit :[font = input; nowordwrap; ] beleg = Apply[Plus, x] :[font = text; inactive; ] Zur Probe wird noch einmal der Stichprobenumfang berechnet. :[font = input; nowordwrap; ] spumfang = Sum[beleg[[i]], {i, Length[beleg]}] :[font = text; inactive; ] Die Liste der Ursprungswerte enthält noch jeweils in der ersten Spalte die Zahl der Nichtraucher und Nichtraucherinnen, die für die folgenden Betrachtungen nicht mehr relevant sind. Mit :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] Drop[ Liste, n ] :[font = text; inactive; ] werden die ersten n-Elemente aus der Liste entfernt. Damit ergibt sich der Ursprungsdatensatz für die rauchende Belegschaft mit :[font = input; nowordwrap; ] belegRauch = Drop[Apply[Plus, x], 1] :[font = text; inactive; ] Damit ändert sich auch der Länge des nun zu bearbeitenden Datensatzes auf :[font = input; nowordwrap; ] spumfang = Sum[belegRauch[[i]], {i, Length[belegRauch]} ] :[font = text; inactive; ] Für die folgenden Berechnungen wird auf verschiedene Packages aus dem Bereich Statistik zurückgegriffen, die alle unter dem Verzeichnis ....\statistics\ liegen. Durch das Laden des zugehörigen Package Master.m sind alle Funktionen der in diesem Verzeichnis enthaltenen Packages verfügbar. :[font = input; nowordwrap; ] << Statistics`Master` :[font = subsubsection; inactive; startGroup; ] Lageparameter :[font = text; inactive; ] Die Lageparameter Mittelwert und Median werden mit den nun verfügbaren Funktionen :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] Mean[ Daten ] , Median[ Daten ] :[font = text; inactive; ] berechnet. :[font = input; endGroup; nowordwrap; ] {mittel, med} = {Mean[belegRauch], Median[belegRauch]}; TableForm[N[{mittel, med}], TableHeadings -> { {"Mittelwert:","Median:"}, {"Rauchende Belegschaft"}}] :[font = subsubsection; inactive; startGroup; ] Streuungsparameter :[font = text; inactive; ] Für die wichtigsten Streuungsparameter wird in Mathematica ebenfalls eine Funktion angeboten: :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] Dispersion[ Daten ] :[font = text; inactive; ] Damit werden die Einzelwerte Variance[ ] StandardDeviation[ ] SampleRange[ ] MeanDeviation[ ] MedianDeviation[ ] QuartileDeviation[ ] in einem Aufruf berechnet. :[font = input; nowordwrap; ] streupar = DispersionReport[beleg]; N[streupar] :[font = text; inactive; ] Den nicht implementierten Variationskoeffizienten :[font = special2; inactive; nowordwrap; nohscroll; center; ] v = s / x_bar :[font = text; inactive; ] berechnet man mit :[font = input; endGroup; endGroup; endGroup; nowordwrap; ] varkoef = (Variance /. streupar) / (Mean[beleg]); N[varkoef] :[font = section; inactive; startGroup; ] Lineare Regression :[font = subsubsection; inactive; startGroup; Cclosed; ] Ausgangssituation :[font = text; inactive; ] Der Firmenleitung ist aus verschiedenen Publikationen bekannt, daß Lärm als ein Streß verursachender Faktor einen Einfluß auf das Rauch- verhalten haben kann. Sie beauftragt daher die statistische Abteilung, diesen Einflußfaktor quantitativ zu bewerten. Die statistische Abteilung erhebt daraufhin bei 50 repräsentativ ausgesuchten Beschäftigten (aus dem Produktionsbereich) folgende Daten: :[font = text; inactive; endGroup; ] - täglicher Zigarettenkonsum am Arbeitsplatz [Stück/Tag] - durchschnittlicher täglicher Lärmpegel am Arbeitsplatz [dB] :[font = subsubsection; inactive; startGroup; Cclosed; ] Regressionsanalyse :[font = text; inactive; ] Um einen ersten Überblick über die Art des funktionalen Zusammenhangs zwischen den Variablen Lärm und Zigarettenkonsum zu erhalten, fertigt man einen Scatterplot an. :[font = input; nowordwrap; ] regdat = ReadList["reg.dat", Number, RecordLists -> True] :[font = input; nowordwrap; ] reglistplot = ListPlot[regdat, Prolog -> AbsolutePointSize[2], PlotRange -> {{55, 90}, {0,20} }, AxesLabel -> {"Laerm[dB]","Zig./Tag"}, PlotLabel -> "Scatterplot" ] :[font = text; inactive; ] Ein nach dem Scatterplot zu vermutender linearer Zusammenhang zwischen dem Lärm (erklärende Variable, Regressor) und dem Zigaretten- konsum (zu erklärende Variabe) wird in folgendem Ansatz modelliert: :[font = special2; inactive; nowordwrap; nohscroll; center; ] y_i = beta0 + beta1 x1 + eps_i mit e_i ~ N(0,1) :[font = text; inactive; ] Die Funktionen zur Regressionsanalyse sind in Mathematica in dem Package Statistics`LinearRegression` enthalten. Diese Funktionen stehen durch das Einlesen des Packages Master.m aus dem Verzeichnis ...\statistics\ schon zur Verfügung. Mit :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] Regress[ Daten, {1, Funktionsform}, x ] :[font = text; inactive; ] wird ein Modell der Form Funktionsform mit der Variablen x den Daten angepaßt. Für die Datenmatrix Daten wird angenommen, daß die letzte Spalte (das letzte Element jeder Liste) die Werte der abhängigen Variablen enthält. :[font = text; inactive; ] Das Symbol x bekam weiter oben schon einen Wert zugewiesen. Diese Zuweisung und x an sich ohne Zuweisung wird mit ClearAll[ ] gelöscht. :[font = input; nowordwrap; ] ClearAll[x] :[font = input; nowordwrap; ] Regress[regdat, {1, x}, x] :[font = text; inactive; ] Mit der Option OutputControl -> NoPrint und eventuell darauffolgende Variablennamen kann die Ausgabe der Funktion Regress[ ] gezielt gesteuert werden. :[font = input; nowordwrap; ] output = Regress[regdat, {1, ex}, ex, OutputControl -> NoPrint, \ OutputList -> {FitResiduals,BestFitCoefficients}]; :[font = text; inactive; ] Die Funktionsform erhält man mit :[font = special1; inactive; output; nowordwrap; nohscroll; center; ] Fit[ Daten, Funktionsform, x ] :[font = input; nowordwrap; ] func[x_] = Fit[regdat, {1, x}, x] :[font = text; inactive; ] Die graphische Darstellung der Regressionsgeraden ergibt eine erst einmal optisch nachvollziehbare Modellwahl. :[font = input; nowordwrap; ] regplot = Plot[func[x], {x, Min[Transpose[regdat] [[1]] ], Max[Transpose[regdat] [[1]] ]}, DisplayFunction -> Identity ] :[font = input; nowordwrap; ] regall = Show[reglistplot, regplot, PlotLabel -> "Beob., Regression", AxesLabel -> {"Laerm[dB]","Zig./Tag"} ] :[font = text; inactive; ] Die Gültigkeit des geschätzten Regressionsmodells hängt im Wesentlichen von der Übereinstimmung der theoretischen Vorgaben (Modellanahmen) mit den vorliegenden Daten ab. Als Beispiel sei dazu die Analyse der Residuen betrachtet. Ein Plot der Residuen übepr