Variance Covariance Matrix from a N*2 Matrix

[Long Hong]

Dear all,

I am writing to ask something about generating a simple variance covariance matrix from a Nx2 matrix. My codes are as follows:

set obs 50
scalar beta0 = 3
scalar beta1 = 0.5
scalar sigma2 = 9
scalar k = 500
local i = 1
while `i' <= k {
  gen epsilons_`i' = rnormal(0, sqrt(sigma2))
  gen yvals_`i' = beta0 + beta1*xvals + epsilons_`i'
  reg yvals_`i' xvals
  matrix beta = (nullmat(beta) \ e(b))
  drop epsilons_`i' yvals_`i'
  local i = `i' + 1
}

After I get the matrix beta (a matrix with 500 repeatedly estimated coefficient matrix), I don't know how to calculate the variance covariance matrix from the matrix. From what I know in R, it's simply -var(beta)-, but I do not know how to make it work in Stata.

I look forward to hearing from you! Thank you very much in advance!

Best,
Long

[Clyde Schechter]

So it looks like you are trying to simulate regression results where the underlying model is in fact a linear model with stipulated parameters and normal residuals. You are capturing the coefficients from each replication in a matrix, but want to get their covariance. There may be a way to do this directly from your matrix beta using some commands in Mata. But here is a way to do it in Stata directly:

capture program drop mymodel
program define mymodel, eclass
    args beta0 beta1 sigma2
    confirm number `beta0'
    confirm number `beta1'
    assert `sigma2' > 0
    replace y = `beta0' + `beta1'*x + rnormal(0, sqrt(`sigma2'))
    regress y x
    exit
end

// INITIALIZE PARAMETERS
local beta0     3
local beta1     0.5
local sigma2   9
local nreps     500

//    CREATE A DATA SHELL
clear
set obs 50
set seed 1234 // OR YOUR LUCKY NUMBER HERE

// INITIALIZE A DATA SET WITH APPROPRIATE X-VALUES
// YOU DON'T SHOW US WHERE YOUR X's COME FROM
// SO YOU"LL NEED TO FILL THAT IN HERE

gen x = runiform() // USE YOUR ACTUAL CODE TO GET THE X VALUES HERE
gen y = . // CREATE AN EMPTY VARIABLE TO HOLD THE SIMULATED Ys

//    TEMPORARY FILE TO HOLD REGRESSION COEFFICIENTS
tempfile betas

//    RUN THE SIMULATION
simulate _b, saving(`betas') reps(`nreps'): mymodel `beta0' `beta1' `sigma2'

//    BRING THE OUTPUT INTO MEMORY
use `betas', clear

mkmat _b*, matrix(BETA) // SAVE COEFFICIENTS IN A MATRIX

matrix accum CV = _b*, noconstant deviations // CALCULATE COVARIANCE MATRIX
matrix list CV

[Long Hong]

Dear Prof. Schechter,

Thank you very much for your detailed reply!
I thought of using mat accum before, but it produced confusing results. I replicated your codes and the covariance matrix is also very confusing.

it is extremely large while the actual number of the data are very small in magnitudes.

Best regards
Long

[Clyde Schechter]

Oh, sorry. -matrix accum- calculates sums of squares and cross-products, not mean squares and cross-products. So it needs a scale factor:

capture program drop mymodel
program define mymodel, eclass
    args beta0 beta1 sigma2
    confirm number `beta0'
    confirm number `beta1'
    assert `sigma2' > 0
    replace y = `beta0' + `beta1'*x + rnormal(0, sqrt(`sigma2'))
    regress y x
    exit
end

// INITIALIZE PARAMETERS
local beta0     3
local beta1     0.5
local sigma2   9
local nreps     500

//    CREATE A DATA SHELL
clear
set obs 50
set seed 1234 // OR YOUR LUCKY NUMBER HERE

// INITIALIZE A DATA SET WITH APPROPRIATE X-VALUES
// YOU DON'T SHOW US WHERE YOUR X's COME FROM
// SO YOU"LL NEED TO FILL THAT IN HERE

gen x = runiform() // USE YOUR ACTUAL CODE TO GET THE X VALUES HERE
gen y = . // CREATE AN EMPTY VARIABLE TO HOLD THE SIMULATED Ys

//    TEMPORARY FILE TO HOLD REGRESSION COEFFICIENTS
tempfile betas

//    RUN THE SIMULATION
simulate _b, saving(`betas') reps(`nreps'): mymodel `beta0' `beta1' `sigma2'

//    BRING THE OUTPUT INTO MEMORY
use `betas', clear

matrix accum CV = _b*, noconstant deviations // CALCULATE COVARIANCE MATRIX
matrix CV = (1/(`nreps'-1)) * CV
matrix list CV

will get you the covariance matrix.

[Long Hong]

Dear Prof. Schechter, Thank you! It returns the correct matrix now! :D