Untitled

Consider the following estimator

\[\hat{\boldsymbol{\theta}} = \operatorname{argmin}_{\boldsymbol{\theta}} \; \left(\mathbf{X}\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}\right)^T \boldsymbol{\Omega} \left(\mathbf{X}\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}\right),\]

subject to

\[\mathbf{X} \hat{\boldsymbol{\theta}} \geq \hat{\boldsymbol{\nu}}.\]

Assuming \(\boldsymbol{\Omega}\) to be a diagonal matrix (with positive entries) we have

\[\left(\mathbf{X}\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}\right)^T \boldsymbol{\Omega} \left(\mathbf{X}\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}\right) = \left(\boldsymbol{\Omega}^{1/2} \mathbf{X}\boldsymbol{\theta} - \boldsymbol{\Omega}^{1/2} \hat{\boldsymbol{\nu}}\right)^T \left(\boldsymbol{\Omega}^{1/2} \mathbf{X}\boldsymbol{\theta} - \boldsymbol{\Omega}^{1/2} \hat{\boldsymbol{\nu}}\right) = \left(\mathbf{X}^*\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}^*\right)^T \left(\mathbf{X}^*\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}^*\right),\] where \(\mathbf{X}^* = \boldsymbol{\Omega}^{1/2} \mathbf{X}\) and \(\hat{\boldsymbol{\nu}}^* = \boldsymbol{\Omega}^{1/2} \hat{\boldsymbol{\nu}}\). Then, using the results of Liew, 1976, we have

where \(\tilde {\boldsymbol{\theta}}\) corresponds to the standard GMWM estimator, i.e.

\[\begin{align} \hat{\boldsymbol{\theta}} &= \tilde {\boldsymbol{\theta}} + \left[\left(\mathbf{X}^*\right)^{T} \mathbf{X}^*\right]^{-1} \mathbf{X}^T \left\{\mathbf{X} \left[\left(\mathbf{X}^*\right)^{T} \mathbf{X}^*\right]^{-1}\mathbf{X}^T\right\}^{-1} \left(\hat{\boldsymbol{\nu}} - \mathbf{X} \hat{\boldsymbol{\theta}}\right)\\ &= \tilde {\boldsymbol{\theta}} + \left(\mathbf{X}^{T} \boldsymbol{\Omega} \mathbf{X}\right)^{-1} \mathbf{X}^T \left\{\mathbf{X} \left(\mathbf{X}^{T} \boldsymbol{\Omega} \mathbf{X}\right)^{-1}\mathbf{X}^T\right\}^{-1} \left(\hat{\boldsymbol{\nu}} - \mathbf{X} \hat{\boldsymbol{\theta}}\right)\\ \end{align}\]

\[\tilde {\boldsymbol{\theta}} = \operatorname{argmin}_{\boldsymbol{\theta}} \; \left(\mathbf{X}\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}\right)^T \boldsymbol{\Omega} \left(\mathbf{X}\boldsymbol{\theta} - \hat{\boldsymbol{\nu}}\right) = \left[\left(\mathbf{X}^*\right)^{T} \mathbf{X}^*\right]^{-1} \left(\mathbf{X}^*\right)^{T} \hat{\boldsymbol{\nu}}^* = \left(\mathbf{X}^{T} \boldsymbol{\Omega} \mathbf{X}\right)^{-1} \mathbf{X}^T \boldsymbol{\Omega} \hat{\boldsymbol{\nu}}\]

The theoretical WV is given by

\[{\nu}_j = \left[ \frac{1}{\tau_j} \;\;\; \frac{\left(\tau_j^2 + 2\right)}{12 \tau_j} \right] \boldsymbol{\theta}.\]

Let’s simulate a process

library(wv)
library(simts)
n = 10^6
theta = c(1, 0.0000001)
model = WN(sigma2 = theta[1]) + RW(gamma2 = theta[2])
set.seed(21453)
Xt = gen_gts(n = n, model = model)
wv_Xt = wvar(Xt)
plot(wv_Xt)

Let’s add the theoretical WV:

J = floor(log2(n)) - 1
tau = 2^(1:J)
X1 = 1/tau
X2 = (tau^2 + 2)/(12*tau)
X = cbind(X1, X2)
wv_theo = X%*%theta
plot(wv_Xt, legend_position = NA)
lines(tau, wv_theo, col = "red", lwd = 2)

Let’s compute the standard GMWM estimator

Omega = diag(1/(wv_Xt$ci_high - wv_Xt$ci_low)^2)
estim_gmwm = solve(t(X)%*%Omega%*%X)%*%t(X)%*%Omega%*%wv_Xt$variance
plot(wv_Xt, legend_position = NA)
lines(tau, X%*%estim_gmwm, col = "green3", lwd = 2)

Let’s compute the conservative GMWM estimator

XtOmX_inv = solve(t(X)%*%Omega%*%X)
plot(wv_Xt$variance - X%*%estim_gmwm)

index = c(1,2)
A = X[index,]
add = matrix(NA, (J-1), 2)
for (i in 2:J){
  index = c(i)
  A = X[index,]
  add[i-1, ] = XtOmX_inv%*%A%*%solve(A%*%XtOmX_inv%*%A)%*%(wv_Xt$variance[index] - A%*%estim_gmwm)
}
estim_cgmwm = estim_gmwm + apply(add, 2, max)



plot(wv_Xt, legend_position = NA)
lines(tau, X%*%estim_cgmwm, col = "orange2", lwd = 2)

plot(X%*%estim_cgmwm - wv_Xt$variance)