Wie geht glmnet mit Überdispersion um?

9

Ich habe eine Frage zum Modellieren von Text über Zähldaten, insbesondere zum Verwenden der lassoTechnik zum Reduzieren von Features.

Angenommen, ich habe N Online-Artikel und die Anzahl der Seitenaufrufe für jeden Artikel. Ich habe 1 Gramm und 2 Gramm für jeden Artikel extrahiert und wollte eine Regression über die 1,2 Gramm durchführen. Da die Merkmale (1,2 Gramm) weit über der Anzahl der Beobachtungen liegen, wäre das Lasso eine gute Methode, um die Anzahl der Merkmale zu verringern. Außerdem habe ich festgestellt, dass glmnetes sehr praktisch ist, eine Lasso-Analyse durchzuführen.

Die Anzahl der Seitenaufrufe ist jedoch überstreut (Varianz> Mittelwert), glmnetbietet jedoch nicht quasipoisson(explizit) oder negative binomialnur poissonfür Zähldaten. Die Lösung, an die ich gedacht habe, besteht darin, log transformdie Zähldaten (eine unter Sozialwissenschaftlern häufig verwendete Methode) zu verwenden und die Antwortvariable grob einer Normalverteilung zu folgen. Als solches könnte ich möglicherweise die Daten mit der Gaußschen Familie unter Verwendung von modellieren glmnet.

Meine Frage lautet also: Ist das angemessen? Oder soll ich nur poisson für glmnetden Fall , glmnetGriffe quasipoisson? Oder gibt es andere R-Pakete, die mit dieser Situation umgehen?

Vielen Dank!

Sonya S.
quelle

Antworten:

14

Kurze Antwort

Überdispersion spielt bei der Schätzung eines Vektors von Regressionskoeffizienten für den bedingten Mittelwert in einem Quasi / Poisson-Modell keine Rolle! Sie werden in Ordnung sein, wenn Sie die Überdispersion hier vergessen, glmnet mit der Poisson-Familie verwenden und sich nur darauf konzentrieren, ob Ihr kreuzvalidierter Vorhersagefehler gering ist.

Die Qualifikation folgt unten.


Poisson, Quasi-Poisson und Schätzfunktionen:

yμxβ

θV.einr(y)=θμθθθββ^ sind für die Quasi- und Poisson-Modelle identisch!

Lassen Sie mich anhand eines Beispiels veranschaulichen (beachten Sie, dass Sie scrollen müssen, um den gesamten Code und die Ausgabe zu sehen):

> library(MASS)
> data(quine) 
> modp <- glm(Days~Age+Sex+Eth+Lrn, data=quine, family="poisson")
> modqp <- glm(Days~Age+Sex+Eth+Lrn, data=quine, family="quasipoisson")
> summary(modp)

Call:
glm(formula = Days ~ Age + Sex + Eth + Lrn, family = "poisson", 
    data = quine)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-6.808  -3.065  -1.119   1.819   9.909  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  2.71538    0.06468  41.980  < 2e-16 ***
AgeF1       -0.33390    0.07009  -4.764 1.90e-06 ***
AgeF2        0.25783    0.06242   4.131 3.62e-05 ***
AgeF3        0.42769    0.06769   6.319 2.64e-10 ***
SexM         0.16160    0.04253   3.799 0.000145 ***
EthN        -0.53360    0.04188 -12.740  < 2e-16 ***
LrnSL        0.34894    0.05204   6.705 2.02e-11 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 2073.5  on 145  degrees of freedom
Residual deviance: 1696.7  on 139  degrees of freedom
AIC: 2299.2

Number of Fisher Scoring iterations: 5

> summary(modqp)

Call:
glm(formula = Days ~ Age + Sex + Eth + Lrn, family = "quasipoisson", 
    data = quine)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-6.808  -3.065  -1.119   1.819   9.909  

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   2.7154     0.2347  11.569  < 2e-16 ***
AgeF1        -0.3339     0.2543  -1.313 0.191413    
AgeF2         0.2578     0.2265   1.138 0.256938    
AgeF3         0.4277     0.2456   1.741 0.083831 .  
SexM          0.1616     0.1543   1.047 0.296914    
EthN         -0.5336     0.1520  -3.511 0.000602 ***
LrnSL         0.3489     0.1888   1.848 0.066760 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for quasipoisson family taken to be 13.16691)

    Null deviance: 2073.5  on 145  degrees of freedom
Residual deviance: 1696.7  on 139  degrees of freedom
AIC: NA

Number of Fisher Scoring iterations: 5

Wie Sie sehen können, deviance(modp)/modp$df.residualändern sich die Regressionskoeffizienten (Punktschätzungen) überhaupt nicht, obwohl wir in diesem Datensatz eine starke Überdispersion von 12,21 haben . Beachten Sie jedoch, wie sich die Standardfehler ändern.

Die Frage nach dem Effekt der Überdispersion in bestraften Poisson-Modellen

θβ

G(μ)=xβ+f(β)

βθμ

glmnetf(β)=0

> library(glmnet)
> y <- quine[,5]
> x <- model.matrix(~Age+Sex+Eth+Lrn,quine)
> modl <- glmnet(y=y,x=x, lambda=c(0.05,0.02,0.01,0.005), family="poisson")
> coefficients(modl)
8 x 4 sparse Matrix of class "dgCMatrix"
                    s0         s1         s2         s3
(Intercept)  2.7320435  2.7221245  2.7188884  2.7172098
(Intercept)  .          .          .          .        
AgeF1       -0.3325689 -0.3335226 -0.3339580 -0.3340520
AgeF2        0.2496120  0.2544253  0.2559408  0.2567880
AgeF3        0.4079635  0.4197509  0.4236024  0.4255759
SexM         0.1530040  0.1581563  0.1598595  0.1607162
EthN        -0.5275619 -0.5311830 -0.5323936 -0.5329969
LrnSL        0.3336885  0.3428815  0.3459650  0.3474745

Was macht OD also mit bestraften Regressionsmodellen? Wie Sie vielleicht wissen, gibt es immer noch einige Debatten darüber, wie Standardfehler für bestrafte Modelle richtig berechnet werden können (siehe z. B. hier ), und sie glmnetwerden wahrscheinlich aus diesem Grund ohnehin nicht ausgegeben. Es kann sehr gut sein, dass die OD den Inferenzteil des Modells beeinflusst, genau wie im nicht bestraften Fall, aber wenn in diesem Fall kein Konsens über die Inferenz erzielt wird, werden wir es nicht wissen.

Abgesehen davon kann man all diese Unordnung hinter sich lassen, wenn man bereit ist, eine Bayes'sche Sichtweise zu vertreten, in der bestrafte Modelle nur Standardmodelle mit einem bestimmten Prior sind.

Momo
quelle
@ Mono, danke für deine sehr ausführliche Erklärung! Hier ist mein Verständnis, und bitte korrigieren Sie mich, wenn ich falsch liege: poissonund quasipoissonRegressionen schätzen Koeffizienten auf die gleiche Weise und unterscheiden sich darin, wie sie Standardfehler und damit die Signifikanz schätzen. Bei der Lasso-Methode ist die Berechnung von Standardfehlern jedoch noch nicht konsensfähig, und daher liegt ihre derzeitige Verwendung hauptsächlich in der Variablenauswahl und nicht in der Inferenz. Als solches spielt es keine Rolle, ob wir glmnetmit Poisson oder Quasipoisson arbeiten, aber was bewirkt, ist, dass der kreuzvalidierte Fehler minimiert werden sollte.
Sonya S.
@Mono, noch eine Anmerkung, ich habe summary(modqp)mich selbst laufen lassen und gesehen, dass es genau die gleichen Koeffizientenschätzungen hat. Ich glaube, Ihre Antwort wird mehr Menschen in dieser Angelegenheit zugute kommen, da ich keine gefunden habe. Ich schlage daher vor, dass Sie die Ausgabe der Zusammenfassung (modqp) hinzufügen, um ein noch besser illustriertes Beispiel zu erhalten. Nochmals vielen Dank!
Sonya S.
1
@Sonya Ihre ist eine gute Zusammenfassung. Der Schlüssel ist, dass bei der Schätzung der Parameter für den bedingten Mittelwert die Schätzfunktionen (z. B. die Score-Funktion) für Poisson und Quasipoisson gleich sind! Für diese Parameter spielt es daher keine Rolle, ob eine Bestrafung vorliegt oder nicht, solange es sich um dieselbe Bestrafung handelt. Ich mache das oben deutlicher. Danke auch für den Zeiger bezüglich der Zusammenfassung (modq), aber der ist schon da, er wird auf einem normalen Bildschirm nur "eingerahmt", also muss man nach unten scrollen.
Momo
Ich frage mich immer noch, ob es möglich ist, dass in Poisson weniger Variablen geschrumpft werden als in einer Quasi-Poisson-Spezifikation, die korrekter ist und wahrscheinlich zu einer besseren Vorhersagegenauigkeit führen würde als das Poisson-Modell, da das Stichprobenmodell korrekter ist.
Dreistes Gleichgewicht
In diesem Sinne könnte es auch sein, dass in Poisson mehr Variablen geschrumpft werden, als in Fällen von UNDER-Dispersion geschrumpft werden sollten (z. B. wenn Sie ein robustes Poisson-Modell verwenden, um die relativen Risikoverhältnisse für 0/1-Daten zu schätzen).
Dreistes Gleichgewicht