Bootstrap-Filter / Partikelfilter-Algorithmus (Verständnis)

12

Ich habe wirklich ein Unverständnis darüber, wie der Bootstrap-Filter funktioniert. Ich kenne die Konzepte grob, kann aber bestimmte Details nicht erfassen. Diese Frage ist für mich, um das Durcheinander zu beseitigen. Hier werde ich diesen beliebten Filteralgorithmus aus einer Referenz von Doucet verwenden (bisher denke ich, dass dies die einfachste Referenz ist). Lassen Sie mich zunächst sagen, dass mein Problem darin besteht, zu verstehen, welche Distributionen bekannt und welche unbekannt sind.

Bootstrap-Filter

Das sind meine Fragen:

  1. Wie ist in 2) die Verteilung ? Ist diese Distribution bekannt ? Kennen wir diese Verteilung für alle t ? Wenn ja, aber was ist, wenn wir nicht davon probieren können? Es ist lustig, dass sie diesen wichtigen Stichprobenschritt nennen, aber ich sehe keine Angebotsverteilung.p(xt|xt1(i))t
  2. Auch in 2) ist eine bekannte Verteilung? „Normalisieren Wichtigkeitsgewichte Mittel w ( i ) t = ~ w ( i ) tp(yt|x~t(i)) ? Was bedeutet die Tilde aufxundw? Bedeutet es so etwas wie nicht erneut abgetastet oder nicht normalisiert?wt(i)=w~t(i)i=1Nw~t(i)xw
  3. Ich würde mich freuen, wenn jemand ein einfaches Spielzeugbeispiel mit bekannten Distributionen zur Verwendung dieses Bootstrap-Filters nennen könnte. Das Endziel des Bootstrap-Filters ist mir nicht klar.
tintinthong
quelle

Antworten:

10
  1. Das ist die Übergangsdichte des Zustands ( ), der Teil Ihres Modells ist und daher bekannt ist. Sie müssen im Basisalgorithmus davon abtasten, aber Annäherungen sind möglich. p ( x t | x t - 1 ) ist in diesem Fall die Angebotsverteilung. Es wird verwendet, weil die Verteilung p ( x t | x 0 : t - 1 , y 1 : t ) im Allgemeinen nicht nachvollziehbar ist.xtp(xt|xt1) p(xt|x0:t1,y1:t)

  2. x~xw~wx

  3. p(xt|y1:t)tt

Betrachten Sie das einfache Modell:

Xt=Xt1+ηt,ηtN(0,1)
X0N(0,1)
Yt=Xt+εt,εtN(0,1)

YXp(Xt|Y1,...,Yt)

Xt|Xt1N(Xt1,1)
X0N(0,1)
Yt|XtN(Xt,1)

Anwendung des Algorithmus:

  1. NX0(i)N(0,1)

  2. X1(i)|X0(i)N(X0(i),1)N

    w~t(i)=ϕ(yt;xt(i),1)ϕ(x;μ,σ2)μσ2yt

  3. wtxx0:t(i)

Kehren Sie zu Schritt 2 zurück und fahren Sie mit der erneut abgetasteten Version der Partikel fort, bis wir die gesamte Serie verarbeitet haben.

Eine Implementierung in R folgt:

# Simulate some fake data
set.seed(123)

tau <- 100
x <- cumsum(rnorm(tau))
y <- x + rnorm(tau)

# Begin particle filter
N <- 1000
x.pf <- matrix(rep(NA,(tau+1)*N),nrow=tau+1)

# 1. Initialize
x.pf[1, ] <- rnorm(N)
m <- rep(NA,tau)
for (t in 2:(tau+1)) {
  # 2. Importance sampling step
  x.pf[t, ] <- x.pf[t-1,] + rnorm(N)

  #Likelihood
  w.tilde <- dnorm(y[t-1], mean=x.pf[t, ])

  #Normalize
  w <- w.tilde/sum(w.tilde)

  # NOTE: This step isn't part of your description of the algorithm, but I'm going to compute the mean
  # of the particle distribution here to compare with the Kalman filter later. Note that this is done BEFORE resampling
  m[t-1] <- sum(w*x.pf[t,])

  # 3. Resampling step
  s <- sample(1:N, size=N, replace=TRUE, prob=w)

  # Note: resample WHOLE path, not just x.pf[t, ]
  x.pf <- x.pf[, s]
}

plot(x)
lines(m,col="red")

# Let's do the Kalman filter to compare
library(dlm)
lines(dropFirst(dlmFilter(y, dlmModPoly(order=1))$m), col="blue")

legend("topleft", legend = c("Actual x", "Particle filter (mean)", "Kalman filter"), col=c("black","red","blue"), lwd=1)

Das resultierende Diagramm:enter image description here

Ein nützliches Tutorial ist das von Doucet und Johansen, siehe hier .

Chris Haug
quelle
X1(i)|X0(i)N(0,1)X1(i)|X0(i)N(X0(i),1)
Tintinthong
Das ist richtig, ich habe den Tippfehler behoben
Chris Haug
Die Pfade müssen nicht neu abgetastet werden, oder? Aus anderer Literatur ist es nicht erforderlich, die Pfade abzutasten. Ich muss nur die Partikel bei jedem Zeitschritt probieren. Ich habe mich gefragt, ob es einen Grund gibt, die Pfade neu
abzutasten