sEM | R Documentation |
Alteranting between expectation step (by kalman filter and smoothing) and maximization step (by generalized Dantzig selector for transiton matrix) to estimate transtion matrix and error variances.
sEM(
Y,
A_init,
sig2_eta_init,
sig2_epsilon_init,
Ti_train = NULL,
Ti_gap = NULL,
tol_seq = NULL,
ht_seq = 0,
is_cv = TRUE,
thres = 0.001,
count_vanish = 1,
n_em = NULL,
is_echo = FALSE,
is_sparse = TRUE
)
Y |
observations of time series, a p by T matrix. |
A_init |
a p by p matrix as initial value of transition matrix |
sig2_eta_init |
scalar; initial value of propagation error variance |
sig2_epsilon_init |
scalar; initial value of measurement error variance |
Ti_train |
scalar; the number of time points in training data in cross-validation. |
Ti_gap |
scalar; the number of time points between test data and train data in cross-validation. |
tol_seq |
vector; grid of tolerance parameter in Dantzig selector for cross-validation. If |
ht_seq |
vector; grid of hard-thresholding levels for transition matrix estimate. If |
is_cv |
logical; if true, run cross-validation to tune Dantzig selector tolerance parameter each sparse EM iteration. |
thres |
scalar; if the difference between updates of two consecutive iterations is less that |
count_vanish |
scalar; if the difference between updates of two consecutive
iterations is less that |
n_em |
scalar; the maximal allowed number of EM iterations, otherwise the algorithm is terminated due to too many iterations.
If |
is_echo |
logical; if true, display the information of CV-optimal (tol, ht) each iteration, and of algorithm termination. |
is_sparse |
logical; if false, use standard EM algorithm, and arguments for cross-validation are not needed. |
a list of parameter estimates.
A_est | estimate of transition matrix A . |
sig2_eta_hat | estimate of propagation error variance \sigma_\eta^2 . |
sig2_epsilon_hat | estimate of measurement error variance \sigma_\epsilon^2 . |
iter_err | the difference between updates of two consecutive iterations. |
iter_err_ratio | the difference ratio (over the previous estimate) between updates of two consecutive iterations. |
Xiang Lyu, Jian Kang, Lexin Li
p= 3; Ti=20 # dimension and time
A=diag(1,p) # transition matrix
sig_eta=sig_epsilon=0.2 # error std
Y=array(0,dim=c(p,Ti)) #observation t=1, ...., Ti
X=array(0,dim=c(p,Ti)) #latent t=1, ...., T
Ti_burnin=30 # time for burn-in to stationarity
for (t in 1:(Ti+Ti_burnin)) {
if (t==1){
x1=rnorm(p)
} else if (t<=Ti_burnin) { # burn in
x1=A%*%x1+rnorm(p,mean=0,sd=sig_eta)
} else if (t==(Ti_burnin+1)){ # time series used for learning
X[,t-Ti_burnin]=x1
Y[,t-Ti_burnin]=X[,t-Ti_burnin]+rnorm(p,mean=0,sd=sig_epsilon)
} else {
X[,t- Ti_burnin]=A%*%X[,t-1- Ti_burnin]+rnorm(p,mean=0,sd=sig_eta)
Y[,t- Ti_burnin]=X[,t- Ti_burnin]+rnorm(p,mean=0,sd=sig_epsilon)
}
}
sEM_fit=sEM(Y,diag(0.5,p),0.1,0.1,Ti*0.5,Ti*0.2,c(0.01,0.1))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.