hdiVAR (version 1.0.2)

Mstep: maximization step of sparse expectation-maximization algorithm for updating error standard deviations

Description

Update \(\sigma_\eta,\sigma_\epsilon\) based on estimate of A and conditional expecation and covariance from expectation step.

Usage

Mstep(Y,A,EXtT,EXtt,EXtt1,is_MLE=FALSE)

Value

a list of estimates of error standard deviations.

sig_eta	estimate of \(\sigma_\eta\).
sig_epsilon	estimate of \(\sigma_\epsilon\).
A	naive MLE of transition matrix \(A\) by conditional expecation and covariance from expecation step. Exist if is_MLE=TRUE.

Arguments

Y

observations of time series, a p by T matrix.

A

current estimate of transition matrix \(A\). If is_MLE=TRUE, use naive MLE of transition matrix, by conditional expecation and covariance from expecation step, to update error standard deviations.

EXtT

a p by T matrix of column \(E[x_t | y_1,\ldots,y_T, \hat{A},\hat{\sigma}_\eta,\hat{\sigma}_\epsilon]\) from expectation step.

EXtt

a p by p by T tensor of first-two-mode slice \(E[x_t x_t^\top | y_1,\ldots,y_T, \hat{A},\hat{\sigma}_\eta,\hat{\sigma}_\epsilon]\) from expectation step.

EXtt1

a p by p by T-1 matrix of first-two-mode slice \(E[x_t x_{t+1}^\top | y_1,\ldots,y_T, \hat{A},\hat{\sigma}_\eta,\hat{\sigma}_\epsilon]\) from expectation step.

is_MLE

logic; if true, use naive MLE of transition matrix, by conditional expecation and covariance from expecation step, to update error variances. Otherwise, use input argument A.

Author

Xiang Lyu, Jian Kang, Lexin Li

Examples

p= 2; Ti=10  # 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=100 # 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)
  }
}

# expectation step
Efit=Estep(Y,A,sig_eta,sig_epsilon,x1,diag(1,p))
EXtT=Efit[["EXtT"]]
EXtt=Efit[["EXtt"]]
EXtt1=Efit[["EXtt1"]]
# maximization step for error standard deviations
Mfit=Mstep(Y,A,EXtT,EXtt,EXtt1)