Regularisierte bayesianische logistische Regression in JAGS

Da die L1-Regularisierung einem Laplace (doppelt exponentiell) vor den relevanten Koeffizienten entspricht, können Sie dies wie folgt tun. Hier habe ich drei unabhängige Variablen x1, x2 und x3 und y ist die binäre Zielvariable. Auswahl des Regularisierungsparameters $\lambda$ erfolgt hier durch Setzen eines Hyperpriors, in diesem Fall gerade gleichmäßig über einen großen Bereich.

model {
  # Likelihood
  for (i in 1:N) {
    y[i] ~ dbern(p[i])

    logit(p[i]) <- b0 + b[1]*x1[i] + b[2]*x2[i] + b[3]*x3[i]
  }

  # Prior on constant term
  b0 ~ dnorm(0,0.1)

  # L1 regularization == a Laplace (double exponential) prior 
  for (j in 1:3) {
    b[j] ~ ddexp(0, lambda)  
  }

  lambda ~ dunif(0.001,10)
  # Alternatively, specify lambda via lambda <- 1 or some such
}

Probieren wir es mit dem dclone Paket in R aus!

library(dclone)

x1 <- rnorm(100)
x2 <- rnorm(100)
x3 <- rnorm(100)

prob <- exp(x1+x2+x3) / (1+exp(x1+x2+x3))
y <- rbinom(100, 1, prob)

data.list <- list(
  y = y,
  x1 = x1, x2 = x2, x3 = x3,
  N = length(y)
)

params = c("b0", "b", "lambda")

temp <- jags.fit(data.list, 
                 params=params, 
                 model="modela.jags",
                 n.chains=3, 
                 n.adapt=1000, 
                 n.update=1000, 
                 thin=10, 
                 n.iter=10000)

Und hier sind die Ergebnisse im Vergleich zu einer unregelmäßigen logistischen Regression:

> summary(temp)

<< blah, blah, blah >> 

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

          Mean     SD Naive SE Time-series SE
b[1]   1.21064 0.3279 0.005987       0.005641
b[2]   0.64730 0.3192 0.005827       0.006014
b[3]   1.25340 0.3217 0.005873       0.006357
b0     0.03313 0.2497 0.004558       0.005580
lambda 1.34334 0.7851 0.014333       0.014999

2. Quantiles for each variable: << deleted to save space >>

> summary(glm(y~x1+x2+x3, family="binomial"))

  << blah, blah, blah >>

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  0.02784    0.25832   0.108   0.9142    
x1           1.34955    0.32845   4.109 3.98e-05 ***
x2           0.78031    0.32191   2.424   0.0154 *  
x3           1.39065    0.32863   4.232 2.32e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

<< more stuff deleted to save space >>

Und wir können sehen, dass die drei bParameter tatsächlich gegen Null geschrumpft sind.

Ich weiß leider nicht viel über Prioritäten für den Hyperparameter der Laplace-Verteilung / den Regularisierungsparameter. Ich neige dazu, gleichmäßige Verteilungen zu verwenden und schaue auf den hinteren Teil, um zu sehen, ob er einigermaßen gut aussieht, z. Bisher war dies normalerweise der Fall. Es funktioniert auch für mich, es als Varianzparameter zu behandeln und die Empfehlungen von Gelman Prior-Verteilungen für Varianzparameter in hierarchischen Modellen zu verwenden.

Bogenschütze
quelle

Du bist der beste! Ich werde die Frage für eine Weile offen lassen, falls jemand eine andere Implementierung hat. Zum einen scheint es, dass binäre Indikatoren verwendet werden können, um variable Ein- / Ausschlüsse zu erzwingen. Dies kompensiert die Tatsache, dass unter Bayesian Lasso die Variablenauswahl nicht tatsächlich erfolgt, da die Betas mit dem doppelten Exponentialprior keine Posterioren haben, die genau Null sind.

Jack Tanner

Richtig, das mache ich auch. Es ist ähnlich wie beim Erstellen eines Priores mit einer Punktmasse von 0 und anschließender Verwendung des Nullentricks zum Abtasten des Priores

b_{i}

$b_i$ wird dann eine Mischung aus einer Punktmasse bei 0 und einem Dublettenexponential, obwohl der Code unterschiedlich ist. Ich bin gespannt, was andere Leute tun, +1 auf die Frage.

Bogenschütze

Regularisierte bayesianische logistische Regression in JAGS

Antworten: