# Calculating variance / standard errors for a Weighted Repeat Sales model

I'm writing an implementation of the Case-Shiller Real Estate Index, which is based on a variation of the weighted least squares, except for the introduction of a dummy matrix Z. I've calculated the index successfully, but I would like to show a 95% confidence interval along with the value. No where in the literature could I find the equation for how to calculate a confidence interval.

The math behind the C-S index is called a Weighted Repeat Sales (WRS) model and is done in three steps (page 19 of [1]):

[First], in order to estimate consistent estimates of the model coefficients, $$\beta$$, we use an instrumental variables estimator, $$\beta = (Z' X )^{−1} Z'Y$$, where Z is a matrix with N rows and T−1 columns that indicates when the sales for each property occurred. The Z matrix is constructed by replacing the positive or negative price levels in X with 1 or –1, respectively.

Then, step two uses the predicted values from this regression to calculate a weights diagonal matrix, which is called $$\Omega$$ and is used in step 3:

If the errors of the value-weighted arithmetic repeat sales model have this heteroskedastic variance structure, then more precise index estimates can be produced by estimating a weighted regression model, $$\beta = (Z' \Omega^{-1} X )^{-1} Z' \Omega^{-1} Y$$, where $$\Omega$$ is a diagonal matrix containing the combined mispricing and interval error variance for each sale pair.

My question is what formula do I use to calculate the standard errors of the $$\beta$$ in this weighted model? I tried calculating variance as $$\hat\sigma^{2} = \frac{Y' \Omega^{-1} (I - H) Y }{n - p}$$ where H is the hat matrix $$H = X (Z' \Omega^{-1} X)^{-1} X' \Omega^{-1}$$ but here I have this new Z variable that should probably be used somewhere, and this formula takes forever to compute and I don't believe produces the right results.

Here is my Julia code in case it's helpful.

Y = -X[:, 1]
X = X[:, 2:end]
Z = sign.(X)

β = inv(Z'X)Z'Y

e = (Y - X * β)[:]
h = inv(Δ' * Δ) * Δ' * e.^2
h = max.(h, 0)

w = sqrt.(Δ * h)[:]
ω = 1 ./ w
Ω = Diagonal(w)

β = inv(Z' * inv(Ω) * X) * Z' * inv(Ω) * Y

# compute confidence interval
H = X * inv(Z' * inv(Ω) * X) * X' * inv(Ω)
σ̂² = Y' * inv(Ω) * (I - H) * Y / (n - p)
Var = σ̂² * inv(X' * inv(Ω) * X)
t = quantile(TDist(p), 0.95)
se = sqrt.(diag(Var))
insert!(β, 1, 1)    # index starts at 1
insert!(se, 1, 0)   # std error of first period is 0
return DataFrame(index=1 ./ β,
ci₀ = 1 ./ (β - t * se),
ci₁ = 1 ./ (β + t * se),
dates=rng)


$$var(\hat\beta- \beta) = (Z'X)^{-1}V(X'Z)^{-1}$$
where $$V=\Sigma{Z'\hat u \hat u' Z}$$ and $$\hat u = Y - X\beta$$.