capture program drop reg_matrix
program define reg_matrix, eclass

  version 8.0
 
    syntax varlist(min=2 numeric) [if] [in] [, Level(integer $S_level)]
    marksample touse                       /* mark cases in the sample */
    tokenize "`varlist'"

    quietly matrix accum sscp = `varlist' if `touse'
    local nobs = r(N)
    local df = `nobs' - (rowsof(sscp) - 1) /* df residual */

    matrix XX = sscp[2...,2...]            /* X'X */
    matrix Xy = sscp[1,2...]               /* X'y */

    matrix b = Xy * syminv(XX)             /* (X'X)-1X'y */
    local k = colsof(b)                    /* number of coefs */
    matrix hat = Xy * b'
    matrix V = syminv(XX) * (sscp[1,1] - hat[1,1])/`df'
    matrix C = corr(V)
    matrix seb = vecdiag(V)
    matrix seb = seb[1, 1...]
    matrix t = J(1,`k',0)
    matrix p = t

    local i = 1
    while `i' <= `k' {
      matrix seb[1,`i'] = sqrt(seb[1,`i'])
      matrix t[1,`i'] = b[1,`i']/seb[1,`i']
      matrix p[1,`i'] = tprob(`df',t[1,`i'])
      local i = `i' + 1
    }

 ereturn post b V, dof(`df') obs(`nobs') depname(`1') /*
      */ esample(`touse')
    ereturn local depvar "`1'"
    ereturn local cmd "matreg"
    
    display
    ereturn display, level(`level')
    
  display
  display "Covariance of the regression coefficients"
  matrix list e(V)
  display
  display "Correlation of the regression coefficients"
  matrix list C

  matrix drop sscp XX Xy hat seb t p C 
  ereturn clear
end