Newey-West standard errors in OLS I am trying to compute robust coefficient estimates for OLS, using the hac() function in MATLAB (see description of function in MathWorks).
In my case, I am regressing a Y variable against a constant only (i.e. X is a vector of ones). In my real example, the data is highly autocorrelated, hence the importance of doing HAC adjustment. Now, here is a simple example (note that Y here is not autocorrelated):
Y = rand(500,1);
X = ones(500, 1);  
hac(X, Y, 'intercept', false, 'weights','BT','display','full')

However, when I compare the results to simple OLS regression:
lscov(X, Y)

it turns out that the coefficient (here intercept) tandard error is exactly the same, as if no HAC adjustment was done when using hac(). Am I doing something wrong here?
Any help would be appreciated! Thank you!
 A: MATLAB's hac() function adjusts your variance-covariance matrix to be consistent in the presence of heteroskedasticity or autocorrelation. The coefficients do not change because OLS is still consistent in the presence of heteroskedasticity since it does not use any covariances to estimate parameters. This is the first line from MATLAB help on the function

Heteroscedasticity and autocorrelation consistent covariance
  estimators

UPDATE
The problem is in your code. You're not testing hac() properly. Its purpose is to adjust the covariances of parameters to the correlation in errors. So, to test this you have to generate correlated errors. Look at this code, and compare with yours:
r = rand(500,1);
e = zeros(500,1);
phi = 0.5;
e(2:end) = phi*r(2:end) + (1-phi)*r(1:end-1);

%%

X = ones(500, 1);  
Y = X + e;
hac(X, Y, 'intercept', false, 'weights','BT','display','full')

You'll see the difference in results immediately.
A: That would indeed be surprising, and although @Aksakal is right that without serial correlation, the estimates should not differ too much, it would still be a fluke if they were exactly equal (and I presume you have run your code more than once). 
One possible explanation (hence at best a partial answer) could be that your bandwidth (I am not sure how they are precisely chosen) is such that your HAC estimates boil down to HC estimates, which, in turn, boil down to the classical ones when only regressing on a constant:
The robust variance matrix is
$\hat V= n(X'X)^{-1}(X'\hat\Omega X)(X'X)^{-1}$
with $\hat\Omega$ the covariance matrix of the residuals. Here, $X$ is just a column of 1s. Hence, $(X'X)^{-1}=1/n$. Now, if $\hat\Omega$ is a diagonal matrix (possibly because the weights are such that only the off-diagonal elements receive zero weight), we obtain
$$X'\hat\Omega X=\sum_{i=1}^n\hat u_i^2$$ 
such that
$\hat V= n\frac{1}{n}\sum_{i=1}^n\hat u_i^2\frac{1}{n}=\frac{1}{n}\sum_{i=1}^n\hat u_i^2$
Upon inserting the degrees of freedom correction $n/(n-1)$ used in the default OLS formula, we get
$$\frac{1}{n-1}\sum_{i=1}^n\hat u_i^2.$$ In R (I do not have access to MATLAB):
library(sandwich)
n <- 10
y <- rnorm(n,1)
reg <- lm(y~1)
vcov(reg)
vcovHC(reg,type="HC0")*n/(n-1)

