Montage d'un coefficient DLM variable dans le temps

9

Je veux adapter un DLM avec des coefficients variant dans le temps, c'est-à-dire une extension de la régression linéaire habituelle,

.yt=θ1+θ2x2

J'ai un prédicteur ( ) et une variable de réponse ( y t ), des captures annuelles de poissons marins et intérieurs respectivement de 1950 à 2011. Je veux que le modèle de régression DLM suive,x2yt

yt=θt,1+θt,2xt

où l'équation d'évolution du système est

θt=Gtθt1

à partir de la page 43 de Modèles linéaires dynamiques avec R par Petris et al.

Du codage ici,

fishdata <- read.csv("http://dl.dropbox.com/s/4w0utkqdhqribl4/fishdata.csv", header=T)
x <- fishdata$marinefao
    y <- fishdata$inlandfao

lmodel <- lm(y ~ x)
summary(lmodel)
plot(x, y)
abline(lmodel)

Il est clair que les coefficients variant dans le temps du modèle de régression sont plus appropriés ici. Je suis son exemple des pages 121 - 125 et je veux l'appliquer à mes propres données. Ceci est le codage de l'exemple

############ PAGE 123
require(dlm)

capm <- read.table("http://shazam.econ.ubc.ca/intro/P.txt", header=T)
capm.ts <- ts(capm, start = c(1978, 1), frequency = 12)
colnames(capm)
plot(capm.ts)
IBM <- capm.ts[, "IBM"]  - capm.ts[, "RKFREE"]
x <- capm.ts[, "MARKET"] - capm.ts[, "RKFREE"]
x
plot(x)
outLM <- lm(IBM ~ x)
outLM$coef
    acf(outLM$res)
qqnorm(outLM$res)
    sig <- var(outLM$res)
sig

mod <- dlmModReg(x,dV = sig, m0 = c(0, 1.5), C0 = diag(c(1e+07, 1)))
outF <- dlmFilter(IBM, mod)
outF$m
    plot(outF$m)
outF$m[ 1 + length(IBM), ]

########## PAGES 124-125
buildCapm <- function(u){
  dlmModReg(x, dV = exp(u[1]), dW = exp(u[2:3]))
}

outMLE <- dlmMLE(IBM, parm = rep(0,3), buildCapm)
exp(outMLE$par)
    outMLE
    outMLE$value
mod <- buildCapm(outMLE$par)
    outS <- dlmSmooth(IBM, mod)
    plot(dropFirst(outS$s))
outS$s

Je veux pouvoir tracer les estimations de lissage plot(dropFirst(outS$s))pour mes propres données, que j'ai du mal à exécuter.

MISE À JOUR

Je peux maintenant produire ces tracés mais je ne pense pas qu'ils soient corrects.

fishdata <- read.csv("http://dl.dropbox.com/s/4w0utkqdhqribl4/fishdata.csv", header=T)
x <- as.numeric(fishdata$marinefao)
    y <- as.numeric(fishdata$inlandfao)
xts <- ts(x, start=c(1950,1), frequency=1)
xts
yts <- ts(y, start=c(1950,1), frequency=1)
yts

lmodel <- lm(yts ~ xts)
#################################################
require(dlm)
    buildCapm <- function(u){
  dlmModReg(xts, dV = exp(u[1]), dW = exp(u[2:3]))
}

outMLE <- dlmMLE(yts, parm = rep(0,3), buildCapm)
exp(outMLE$par)
        outMLE$value
mod <- buildCapm(outMLE$par)
        outS <- dlmSmooth(yts, mod)
        plot(dropFirst(outS$s))

> summary(outS$s); lmodel$coef
       V1              V2       
 Min.   :87.67   Min.   :1.445  
 1st Qu.:87.67   1st Qu.:1.924  
 Median :87.67   Median :3.803  
 Mean   :87.67   Mean   :4.084  
 3rd Qu.:87.67   3rd Qu.:6.244  
 Max.   :87.67   Max.   :7.853  
 (Intercept)          xts 
273858.30308      1.22505 

L'estimation du lissage d'interception (V1) est loin du coefficient de régression lm. Je suppose qu'ils devraient être plus proches les uns des autres.

hgeop
la source

Réponses:

2

Quel est exactement votre problème?

Le seul écueil que j'ai trouvé est que, apparemment,

fishdata <- read.csv("http://dl.dropbox.com/s/4w0utkqdhqribl4,
                     fishdata.csv", header=T)

lit les données sous forme d'entiers. Je devais les convertir pour flotter,

x <- as.numeric(fishdata$marinefao)
y <- as.numeric(fishdata$inlandfao)

avant de pouvoir invoquer les fonctions dlm *.

F. Tusell
la source
Merci pour vos suggestions @F. Tusell; J'ai mis à jour ma question. Les estimations de lissage produites ne sont pas proches des lmodel$coefestimations. Je suppose que les parcelles sont incorrectes mais je peux me tromper.
hgeop
1
Il n'y a aucune raison de s'attendre à ce que les estimations lissées de la pente et de l'ordonnée à l'origine soient proches des bêtas fixes dans la régression linéaire. En particulier, la pente devrait fluctuer énormément.
F. Tusell