komische Dummy-Koeffizientennamen

Modelle zur Korrelations- und Regressionsanalyse

Moderator: EDi

Antworten
bigben
Beiträge: 2771
Registriert: Mi Okt 12, 2016 9:09 am

komische Dummy-Koeffizientennamen

Beitrag von bigben »

Hallo Forum!

Ein Doktorand bei uns hat eine Reihe von Leuten nach ihrem Ausbildungsgrad befragt und danach in Sekunden gemessen, wie lange etwas dauert.
Ich habe also einen ordered factor für die Schulerfahrung und einen numerischen in Sekunden geeicht. Die Daten sehen etwa so aus:

Code: Alles auswählen

messung <- structure(list(edu = structure(c(6L, 4L, 4L, 6L, 2L, 6L, 6L, 
4L, 7L, 9L, 4L, 6L, 4L, 3L, 6L, 8L, 6L, 4L, 6L, 4L, 1L, 6L, 6L, 
6L, 9L, 6L, 6L, 6L, 5L, 6L, 9L, 9L, 3L, 6L, 4L, 4L, 5L, 6L, 6L, 
4L, 4L, 6L, 6L, 6L, 4L, 6L, 6L, 4L, 4L, 6L, 6L, 4L, 4L, 6L, 6L, 
3L, 6L, 6L, 4L, 4L, 6L, 1L, 4L, 4L, 4L, 6L, 4L, 6L, 4L, 4L, 7L, 
6L, 4L, 6L, 4L, 6L, 4L, 4L, 6L, 6L, 4L, 4L, 4L, 9L, 6L, 6L, 6L, 
7L, 9L, 6L, 4L, 6L, 4L, 7L, 6L, 4L, 4L, 3L, 5L, 4L, 4L, 6L, 7L, 
5L, 6L, 6L, 4L, 4L, 6L, 6L, 6L, 8L, 9L, 2L, 6L), .Label = c("5. Klasse", 
"6. Klasse", "7. Klasse", "8. Klasse", "9. Klasse", "10. Klasse", 
"Abitur", "Fachhochschule", "Abitur + Studium"), class = c("ordered", 
"factor")), duration = c(199, 212, 151, 229, 204, 94, 182, 172, 
96, 126, 217, 135, 176, 156, 174, 146, 232, 224, 183, 274, 236, 
238, 119, 185, 192, 239, 138, 125, 105, 131, 197, 190, 182, 140, 
187, 146, 73, 263, 116, 135, 157, 223, 156, 106, 174, 143, 125, 
138, 268, 145, 138, 273, 213, NA, 191, 154, 111, 86, NA, 241, 
215, 407, 116, 276, 455, 190, 200, 182, 325, 104, 136, NA, 244, 
NA, 220, 122, 351, 228, 125, 142, 240, 319, 211, 257, 276, 207, 
122, 104, 76, 181, 157, 227, 162, 168, 84, 167, 160, 448, NA, 
119, 170, 111, 117, 101, 272, 149, 153, 143, 167, 156, 123, NA, 
148, 119, 250)), .Names = c("edu", "duration"), class = "data.frame", row.names = c(NA, 
-115L))
das sieht dann so aus:

Code: Alles auswählen

> str(messung)
'data.frame':	115 obs. of  2 variables:
 $ edu     : Ord.factor w/ 9 levels "5. Klasse"<"6. Klasse"<..: 6 4 4 6 2 6 6 4 7 9 ...
 $ duration: num  199 212 151 229 204 94 182 172 96 126 ...
> levels(messung[,"edu"])
[1] "5. Klasse"        "6. Klasse"        "7. Klasse"        "8. Klasse"       
[5] "9. Klasse"        "10. Klasse"       "Abitur"           "Fachhochschule"  
[9] "Abitur + Studium"
Ultimativ möchte ich gerne in einer multiplen Regression nachweisen, dass der Bildungsabschluss und welche anderen Prädiktoren (hier nicht wiedergegeben) einen Einfluss auf die Dauer `duration` haben. Ein erster Versuch sieht so aus:

Code: Alles auswählen

> summary(lm(duration ~ edu, data=messung))

Call:
lm(formula = duration ~ edu, data = messung)

Residuals:
    Min      1Q  Median      3Q     Max 
-103.32  -44.96  -10.96   32.68  247.68 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  180.545     12.476  14.472  < 2e-16 ***
edu.L       -118.351     42.558  -2.781  0.00648 ** 
edu.Q         91.211     33.945   2.687  0.00844 ** 
edu.C         -6.893     34.591  -0.199  0.84246    
edu^4         33.675     44.942   0.749  0.45544    
edu^5        -57.514     42.866  -1.342  0.18273    
edu^6         70.856     41.257   1.717  0.08900 .  
edu^7        -35.684     27.421  -1.301  0.19612    
edu^8        -56.069     27.116  -2.068  0.04124 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 65.11 on 100 degrees of freedom
  (6 observations deleted due to missingness)
Multiple R-squared:  0.2331,	Adjusted R-squared:  0.1717 
F-statistic: 3.799 on 8 and 100 DF,  p-value: 0.0006249
Kann mir bitte jemand erklären, was die verschiedenene Koeffizientennamen bedeuten und/oder, nach welcher Regel R die bildet :?: Bislang bin ich immer nur über selbsterklärende Dummy-Koeffizientennamen gestolpert, aber normalerweise beginnen die Levelbezeichnungen auch nicht mit Ziffern.

Eine Anschlussfrage stelle ich dann später noch.

LG,
Bernhard
---
Programmiere stets so, dass die Maxime Deines Programmierstils Grundlage allgemeiner Gesetzgebung sein könnte
Benutzeravatar
EDi
Beiträge: 1599
Registriert: Sa Okt 08, 2016 3:39 pm

Re: komische Dummy-Koeffizientennamen

Beitrag von EDi »

Für ordered factors sind die default Kontraste "Orthogonal Polynomial Coding", deshalb die Linear, Quadradic, Cubic, ^4, ...
Bitte immer ein reproduzierbares Minimalbeispiel angeben. Meinungen gehören mir und geben nicht die meines Brötchengebers wieder.

Dieser Beitrag ist lizensiert unter einer CC BY 4.0 Lizenz
Bild.
bigben
Beiträge: 2771
Registriert: Mi Okt 12, 2016 9:09 am

Re: komische Dummy-Koeffizientennamen

Beitrag von bigben »

Danke EDi, das ist die richtige Antwort. Mein erster Versuch, das durch Nachzulesen zu verstehen, ist an meiner Unfähigkeit gescheitert. Gefunden habe ich
This type of coding system should be used only with an ordinal variable in which the levels are equally spaced.
Eine Bedingung, die ich in diesem Fall nicht unterstellen würde, weshalb sich für den aktuellen Fall eine einfache Lösung ergibt, nämlich den geordneten Faktor in einen ungeordneten umzuwandeln und das Nachdenken über coding systems zu verschieben, bis ich mehr Zeit habe, mich damit zu beschäftigen.

Meine Anschlussfrage ist folgende: Eigentlich sind mir das viel zu viele Koeffizienten, denn ich habe nur knapp über hundert Datensätze und möchte möglichst viele Prädiktoren in einer multiplen Regression untersuchen. Es würde mir völlig ausreichen, sagen zu können, dass die Ausbildung eine Rolle spielt (Abstufung, Effektgröße, Richtung - alles egal) und dass ich die Signifikanztests der anderen Prädiktoren mit `edu` korrigiert habe.

Ich dachte mir, dass man das vielleicht mit einem mixed effects modell machen könnte. Wir ergänzen meinen Datensatz oben um weitere Prädiktoren, hier einfach mal um das Geschlecht. Sagen wir

Code: Alles auswählen

messung <- cbind(messung, sex=gl(2,1,115))
Und jetzt zur Prüfung der eigenständigen Signifikanz von sex:

Code: Alles auswählen

library(lme4)
library(lmerTest)
model <- lmer(duration ~ sex + (1|edu), data=messung)
summary(model)
Das ergibt erwartungsgemäß einen Output, in dem sex nicht signifikant ist:

Code: Alles auswählen

Linear mixed model fit by REML t-tests use Satterthwaite approximations to degrees of
  freedom [lmerMod]
Formula: duration ~ sex + (1 | edu)
   Data: messung

REML criterion at convergence: 1218.4

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5668 -0.6762 -0.2937  0.4942  3.7962 

Random effects:
 Groups   Name        Variance Std.Dev.
 edu      (Intercept) 2292     47.87   
 Residual             4328     65.79   
Number of obs: 109, groups:  edu, 9

Fixed effects:
            Estimate Std. Error      df t value Pr(>|t|)    
(Intercept)  179.082     20.458   6.110   8.754 0.000111 ***
sex2           1.834     12.755  98.640   0.144 0.885974    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
     (Intr)
sex2 -0.319
Frage 1: Kann ich das so machen, dass ich zu schätzende Koeffizienten spare, indem ich einen Faktor einsetze wie eine Variable mit random effects? Oder ist das methodisch fragwürdig/unsauber/falsch?

Frage 2: Kann ich dem mixed effects modell irgendwie entnehmen, dass edu tatsächlich relevant ist? Gibt es dafür eine Art Signifikanztest?

LG,
Bernhard
---
Programmiere stets so, dass die Maxime Deines Programmierstils Grundlage allgemeiner Gesetzgebung sein könnte
Benutzeravatar
EDi
Beiträge: 1599
Registriert: Sa Okt 08, 2016 3:39 pm

Re: komische Dummy-Koeffizientennamen

Beitrag von EDi »

Das mit dem ordered factor ist bei genauerem Nachdenken sogar sinnvoll. Dadurch, dass die Levels geordnet sind, kann man einen Trend legen (was bei unordered ja nicht geht). Die polynomialen Kontrasts sind ählich der polynomialen Regression bei kontinuierlichen Variablen und kann man somit eine linear/quadratischen/etc Trend sehen.
Frage 1: Kann ich das so machen, dass ich zu schätzende Koeffizienten spare, indem ich einen Faktor einsetze wie eine Variable mit random effects? Oder ist das methodisch fragwürdig/unsauber/falsch?
Machen kann man das. Ich habe eher die Regel: Fixed Effekts wenn die Levels die ganze Population an Faktorlevels darstellen (gibt also nur die und nur an denen sind wir interessiert). Random Effekts dann wenn die Faktorlevels nur einen Teil (ein sample) der ganzen Population an Faktorlevels ist.
Bei Education findet man für beides Argumente.
Habe vor kurzem einen Vortrag von Andrew Gelman ("Bayesian Data Analysis"...) gesehen, da hatte er solche Ausprägungen/Faktoren (und zwar viele, einschließlich Interaktionen) alles als random intercepts drin. Vielleicht auch eine Unterschied zwischen Fachrichtungen (Gelmen ist ja Sozialwissenschaftler, ich habe mein Wissen eher aus den Naturwissenschaften und da hat man andere Probleme/Experimente/etc).

Meiner Meinung nach kann man das machen (sex würde aber z.b. nicht gehen, da nur zwei Levels).

Zur zweiten Fragen komme ich später nochmal...
Bitte immer ein reproduzierbares Minimalbeispiel angeben. Meinungen gehören mir und geben nicht die meines Brötchengebers wieder.

Dieser Beitrag ist lizensiert unter einer CC BY 4.0 Lizenz
Bild.
Benutzeravatar
EDi
Beiträge: 1599
Registriert: Sa Okt 08, 2016 3:39 pm

Re: komische Dummy-Koeffizientennamen

Beitrag von EDi »

Frage 2: Kann ich dem mixed effects modell irgendwie entnehmen, dass edu tatsächlich relevant ist? Gibt es dafür eine Art Signifikanztest?
Man kann das über simulationen machen. Siehe auch glmm faq.

Code: Alles auswählen

# Models
mod_rand <- lmer(duration ~ sex + (1|edu), data=messung)
mod_0 <- lm(duration ~ sex, data=messung)

# explore model
plot(mod_rand)
fixef(mod_rand)
library(lattice)
dotplot(ranef(mod_rand, condVar = TRUE)) # plot random effects
# no clear trend?


# simulation based Likelihood-Ratio Test for edu
library(RLRsim)
exactLRT(m = pdate(mod_rand, REML = FALSE), # need to refit using ML
         m0 = mod_0)
# see also http://glmm.wikidot.com/faq

Bitte immer ein reproduzierbares Minimalbeispiel angeben. Meinungen gehören mir und geben nicht die meines Brötchengebers wieder.

Dieser Beitrag ist lizensiert unter einer CC BY 4.0 Lizenz
Bild.
Antworten