"""
Diagnostics for SUR and 3SLS estimation
"""
__author__= "Luc Anselin lanselin@gmail.com, \
Pedro V. Amaral pedrovma@gmail.com \
Tony Aburaad taburaad@uchicago.edu"
import numpy as np
import scipy.stats as stats
import numpy.linalg as la
from .sur_utils import sur_dict2mat,sur_mat2dict,sur_corr,spdot
from .regimes import buildR1var,wald_test
__all__ = ['sur_setp','sur_lrtest','sur_lmtest','lam_setp','surLMe','surLMlag']
[docs]def sur_setp(bigB,varb):
'''
Utility to compute standard error, t and p-value
Parameters
----------
bigB : dictionary
of regression coefficient estimates,
one vector by equation
varb : array
variance-covariance matrix of coefficients
Returns
-------
surinfdict : dictionary
with standard error, t-value, and
p-value array, one for each equation
'''
vvb = varb.diagonal()
n_eq = len(bigB.keys())
bigK = np.zeros((n_eq,1),dtype=np.int_)
for r in range(n_eq):
bigK[r] = bigB[r].shape[0]
b = sur_dict2mat(bigB)
se = np.sqrt(vvb)
se.resize(len(se),1)
t = np.divide(b,se)
tp = stats.norm.sf(abs(t))*2
surinf = np.hstack((se,t,tp))
surinfdict = sur_mat2dict(surinf,bigK)
return surinfdict
[docs]def lam_setp(lam,vm):
"""
Standard errors, t-test and p-value for lambda in SUR Error ML
Parameters
----------
lam : array
n_eq x 1 array with ML estimates for spatial error
autoregressive coefficient
vm : array
n_eq x n_eq subset of variance-covariance matrix for
lambda and Sigma in SUR Error ML
(needs to be subset from full vm)
Returns
-------
: tuple
with arrays for standard error, t-value and p-value
(each element in the tuple is an n_eq x 1 array)
"""
vvb = vm.diagonal()
se = np.sqrt(vvb)
se.resize(len(se),1)
t = np.divide(lam,se)
tp = stats.norm.sf(abs(t))*2
return (se,t,tp)
[docs]def sur_lrtest(n,n_eq,ldetS0,ldetS1):
'''
Likelihood Ratio test on off-diagonal elements of Sigma
Parameters
----------
n : int
cross-sectional dimension (number of observations for an equation)
n_eq : int
number of equations
ldetS0 : float
log determinant of Sigma for OLS case
ldetS1 : float
log determinant of Sigma for SUR case (should be iterated)
Returns
-------
(lrtest,M,pvalue) : tuple
with value of test statistic (lrtest),
degrees of freedom (M, as an integer)
p-value
'''
M = n_eq * (n_eq - 1)/2.0
lrtest = n * (ldetS0 - ldetS1)
pvalue = stats.chi2.sf(lrtest,M)
return (lrtest,int(M),pvalue)
[docs]def sur_lmtest(n,n_eq,sig):
'''
Lagrange Multiplier test on off-diagonal elements of Sigma
Parameters
----------
n : int
cross-sectional dimension (number of observations for an equation)
n_eq : int
number of equations
sig : array
inter-equation covariance matrix for null model (OLS)
Returns
-------
(lmtest,M,pvalue) : tuple
with value of test statistic (lmtest),
degrees of freedom (M, as an integer)
p-value
'''
R = sur_corr(sig)
tr = np.trace(np.dot(R.T,R))
M = n_eq * (n_eq - 1)/2.0
lmtest = (n/2.0) * (tr - n_eq)
pvalue = stats.chi2.sf(lmtest,M)
return (lmtest,int(M),pvalue)
[docs]def surLMe(n_eq,WS,bigE,sig):
"""
Lagrange Multiplier test on error spatial autocorrelation in SUR
Parameters
----------
n_eq : int
number of equations
WS : array
spatial weights matrix in sparse form
bigE : array
n x n_eq matrix of residuals by equation
sig : array
cross-equation error covariance matrix
Returns
-------
(LMe,n_eq,pvalue) : tuple
with value of statistic (LMe), degrees
of freedom (n_eq) and p-value
"""
# spatially lagged residuals
WbigE = WS * bigE
# score
EWE = np.dot(bigE.T,WbigE)
sigi = la.inv(sig)
SEWE = sigi * EWE
#score = SEWE.sum(axis=1)
#score.resize(n_eq,1)
# note score is column sum of Sig_i * E'WE, a 1 by n_eq row vector
# previously stored as column
score = SEWE.sum(axis=0)
score.resize(1,n_eq)
# trace terms
WW = WS * WS
trWW = np.sum(WW.diagonal())
WTW = WS.T * WS
trWtW = np.sum(WTW.diagonal())
# denominator
SiS = sigi * sig
Tii = trWW * np.identity(n_eq)
tSiS = trWtW * SiS
denom = Tii + tSiS
idenom = la.inv(denom)
# test statistic
#LMe = np.dot(np.dot(score.T,idenom),score)[0][0]
# score is now row vector
LMe = np.dot(np.dot(score,idenom),score.T)[0][0]
pvalue = stats.chi2.sf(LMe,n_eq)
return (LMe,n_eq,pvalue)
[docs]def surLMlag(n_eq,WS,bigy,bigX,bigE,bigYP,sig,varb):
"""
Lagrange Multiplier test on lag spatial autocorrelation in SUR
Parameters
----------
n_eq : int
number of equations
WS : spatial weights matrix in sparse form
bigy : dictionary
with y values
bigX : dictionary
with X values
bigE : array
n x n_eq matrix of residuals by equation
bigYP : array
n x n_eq matrix of predicted values by equation
sig : array
cross-equation error covariance matrix
varb : array
variance-covariance matrix for b coefficients (inverse of Ibb)
Returns
-------
(LMlag,n_eq,pvalue) : tuple
with value of statistic (LMlag), degrees
of freedom (n_eq) and p-value
"""
# Score
Y = np.hstack((bigy[r]) for r in range(n_eq))
WY = WS * Y
EWY = np.dot(bigE.T,WY)
sigi = la.inv(sig)
SEWE = sigi * EWY
score = SEWE.sum(axis=0) # column sums
score.resize(1,n_eq) # score as a row vector
# I(rho,rho) as partitioned inverse, eq 72
# trace terms
WW = WS * WS
trWW = np.sum(WW.diagonal()) #T1
WTW = WS.T * WS
trWtW = np.sum(WTW.diagonal()) #T2
# I(rho,rho)
SiS = sigi * sig
Tii = trWW * np.identity(n_eq) #T1It
tSiS = trWtW * SiS
firstHalf = Tii + tSiS
WbigYP = WS * bigYP
inner = np.dot(WbigYP.T, WbigYP)
secondHalf= sigi * inner
Ipp= firstHalf + secondHalf #eq. 75
# I(b,b) inverse is varb
# I(b,rho)
bp = sigi[0,] * spdot(bigX[0].T,WbigYP) #initialize
for r in range(1,n_eq):
bpwork = sigi[r,] * spdot(bigX[r].T,WbigYP)
bp = np.vstack((bp,bpwork))
# partitioned part
i_inner = Ipp - np.dot(np.dot(bp.T, varb), bp)
# partitioned inverse of information matrix
Ippi = la.inv(i_inner)
# test statistic
LMlag = np.dot(np.dot(score,Ippi),score.T)[0][0]
# p-value
pvalue = stats.chi2.sf(LMlag,n_eq)
return (LMlag,n_eq,pvalue)
def sur_chow(n_eq,bigK,bSUR,varb):
"""
test on constancy of regression coefficients across equations in
a SUR specification
Note: requires a previous check on constancy of number of coefficients
across equations; no other checks are carried out, so it is possible
that the results are meaningless if the variables are not listed in
the same order in each equation.
Parameters
----------
n_eq : int
number of equations
bigK : array
with the number of variables by equation (includes constant)
bSUR : dictionary
with the SUR regression coefficients by equation
varb : array
the variance-covariance matrix for the SUR regression
coefficients
Returns
-------
test : array
a list with for each coefficient (in order) a tuple with the
value of the test statistic, the degrees of freedom, and the
p-value
"""
kr = bigK[0][0]
test = []
bb = sur_dict2mat(bSUR)
kf = 0
nr = n_eq
df = n_eq - 1
for i in range(kr):
Ri = buildR1var(i,kr,kf,0,nr)
tt,p = wald_test(bb,Ri,np.zeros((df,1)),varb)
test.append((tt,df,p))
return test
def sur_joinrho(n_eq,bigK,bSUR,varb):
"""
Test on joint significance of spatial autoregressive coefficient in SUR
Parameters
----------
n_eq : int
number of equations
bigK : array
n_eq x 1 array with number of variables by equation
(includes constant term, exogenous and endogeneous and
spatial lag)
bSUR : dictionary
with regression coefficients by equation, with
the spatial autoregressive term as last
varb : array
variance-covariance matrix for regression coefficients
Returns
-------
: tuple
with test statistic, degrees of freedom, p-value
"""
bb = sur_dict2mat(bSUR)
R = np.zeros((n_eq,varb.shape[0]))
q = np.zeros((n_eq,1))
kc = -1
for i in range(n_eq):
kc = kc + bigK[i]
R[i,kc] = 1
w,p = wald_test(bb,R,q,varb)
return (w,n_eq,p)