1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
|
#' @export
rarx <- function(x,order=c(1,1,1),lambda=0.95){
y <- outputData(x); u <- inputData(x)
N <- dim(y)[1]
na <- order[1];nb <- order[2]; nk <- order[3]
nb1 <- nb+nk-1 ; n <- max(na,nb1); df <- N-na-nb
yout <- apply(y,2,padZeros,n=n)
uout <- apply(u,2,padZeros,n=n)
uindex <- nk:nb1
if(na!=0) yindex <- 1:na
reg <- function(i) {
# regressor
temp <- numeric(0)
if(na!=0) temp <- c(temp,-yout[i-yindex,])
phi <- c(temp,uout[i-uindex,])
phi
}
# R0 <- reg(n+1)%*%t(reg(n+1))
# Plast <- solve(R0)
Plast <- 10^4*diag(na+nb)
theta <- matrix(0,N,na+nb)
theta[1,] <- Plast%*%reg(n+1)%*%y[1,,drop=FALSE]
yhat <- y
yhat[1,] <- t(reg(n+1))%*%t(theta[1,,drop=FALSE])
for(i in 2:N){
temp <- reg(n+i)
yhat[i,] <- t(temp)%*%t(theta[i-1,,drop=FALSE])
eps_i <- y[i,,drop=FALSE] - yhat[i,,drop=FALSE]
kappa_i <- Plast%*%temp/(lambda+t(temp)%*%Plast%*%temp)[1]
theta[i,] <- t(t(theta[i-1,,drop=F])+eps_i[1]*kappa_i)
Plast <- (diag(na+nb)-kappa_i%*%t(temp))%*%Plast/lambda
}
list(theta=theta,yhat=yhat,P=Plast)
}
|
