Diagnóstico de MCMC Geweke

14

Estoy ejecutando una muestra de Metropolis (C ++) y quiero usar las muestras anteriores para estimar la tasa de convergencia.

Un diagnóstico fácil de implementar que encontré es el diagnóstico de Geweke , que calcula la diferencia entre las dos medias de muestra divididas por su error estándar estimado. El error estándar se estima a partir de la densidad espectral en cero.

Znorte=θ¯UN-θ¯si1norteUNSθUN^(0 0)+1nortesiSθsi^(0 0),

donde , son dos ventanas dentro de la cadena de Markov. Investigué un poco sobre qué son y pero me meto en un lío de literatura sobre densidad espectral de energía y espectro espectral de potencia densidad pero no soy un experto en estos temas; Solo necesito una respuesta rápida: ¿son estas cantidades iguales a la varianza muestral? Si no, ¿cuál es la fórmula para calcularlos?B ^ S A θ ( 0 ) ^ S B θ ( 0 )UNsiSθUN^(0 0)Sθsi^(0 0)

Otra duda sobre este diagnóstico de Geweke es cómo elegir θ ? La literatura anterior dice que es algo funcional θ(X) y debería implicar la existencia de una densidad espectral SθUN^(0 0) , pero por conveniencia supongo que la forma más simple es usar el función de identidad (usar muestras ellos mismos). ¿Es esto correcto?

El paquete R coda tiene una descripción, pero tampoco especifica cómo calcular los valores de S

Yang
fuente
podrías mirar en las entrañas de la codafunción geweke.diagpara ver qué hace ...
Ben Bolker

Respuestas:

4

Puede buscar la geweke.diagfunción en el codapaquete a través del código para ver cómo se calcula la varianza, a través de la llamada a la spectrum.ar0función.


Aquí hay una breve motivación del cálculo de la densidad espectral de un proceso AR ( ) en cero.pag

La densidad espectral de un proceso AR ( ) a frecuenciaλpagλ viene dada por la expresión: donde son los parámetros autorregresivos.

F(λ)=σ2(1-j=1pagαjExp(-2πιjλ))2
αj

Esta expresión se simplifica considerablemente cuando se calcula la densidad espectral de un proceso AR ( ) en : pag0 0

F(0 0)=σ2(1-j=1pagαj)2

El cálculo entonces se vería así (sustituyendo los estimadores habituales por parámetros):

tsAR2 = arima.sim(list(ar = c(0.01, 0.03)), n = 1000)  # simulate an AR(2) process
ar2 = ar(tsAR2, aic = TRUE)  # estimate it with AIC complexity selection

# manual estimate of spectral density at zero
sdMan = ar2$var.pred/(1-sum(ar2$ar))^2

# coda computation of spectral density at zer0
sdCoda = coda::spectrum0.ar(tsAr2)$spec

# assert equality
all.equal(sdCoda, sdMan)
tchakravarty
fuente
0

Consulta la página de wikipedia . Verá , que es la densidad espectral. En su caso, debe usar .SXX(ω)SXX(0 0)

xuhdev
fuente