Merge pull request #138 from pedrovma/1.4.3

Minor adjustments to printouts and spatial impacts

spreg/dgp.py

@@ -291,7 +291,7 @@ def make_wxg(wx,gamma):
 
 def dgp_errproc(u,w,lam=0.5,model='sar',imethod='power_exp'):
  """
- dgp_sar: generates pure spatial error process
+ dgp_errproc: generates pure spatial error process
 
  Arguments:
  ----------
@@ -303,7 +303,7 @@ def dgp_errproc(u,w,lam=0.5,model='sar',imethod='power_exp'):
 
  Returns:
  --------
- y : vector of observations following a pure SAR process
+ y : vector of observations following a spatial AR or MA process
 
  Examples
  --------

spreg/ml_lag.py

@@ -340,9 +340,11 @@ class ML_Lag(BaseML_Lag):
  tolerance criterion in mimimize_scalar function and inverse_product
  spat_diag : boolean
  If True, then compute Common Factor Hypothesis test when applicable
- spat_impacts : boolean
- If True, include average direct impact (ADI), average indirect impact (AII),
- and average total impact (ATI) in summary results
+ spat_impacts : string
+ Include average direct impact (ADI), average indirect impact (AII),
+ and average total impact (ATI) in summary results.
+ Options are 'simple', 'full', 'power', or None.
+ See sputils.spmultiplier for more information.
  vm : boolean
  if True, include variance-covariance matrix in summary
  results
@@ -594,7 +596,7 @@ def __init__(
  slx_lags=0,
  method="full",
  epsilon=0.0000001,
- spat_impacts=True,
+ spat_impacts="simple",
  vm=False,
  spat_diag=True,
  name_y=None,
@@ -642,7 +644,7 @@ def __init__(
  if spat_diag and slx_lags==1:
  diag_out = _spat_diag_out(self, w, 'yend', ml=True)
  if spat_impacts and slx_lags == 0:
- impacts = _summary_impacts(self, _spmultiplier(w, self.rho))
+ impacts = _summary_impacts(self, _spmultiplier(w, self.rho, method=spat_impacts), spat_impacts)
  try:
  diag_out += impacts
  except TypeError:

spreg/ml_lag_regimes.py

@@ -66,9 +66,11 @@ class ML_Lag_Regimes(BaseML_Lag, REGI.Regimes_Frame):
  the spatial parameter is fixed accross regimes.
  spat_diag : boolean
  If True, then compute Common Factor Hypothesis test when applicable
- spat_impacts : boolean
- If True, include average direct impact (ADI), average indirect impact (AII),
- and average total impact (ATI) in summary results
+ spat_impacts : string
+ Include average direct impact (ADI), average indirect impact (AII),
+ and average total impact (ATI) in summary results.
+ Options are 'simple', 'full', 'power', or None.
+ See sputils.spmultiplier for more information.
  cores : boolean
  Specifies if multiprocessing is to be used
  Default: no multiprocessing, cores = False
@@ -318,7 +320,7 @@ def __init__(
  regime_lag_sep=False,
  cores=False,
  spat_diag=True,
- spat_impacts=True,
+ spat_impacts="simple",
  vm=False,
  name_y=None,
  name_x=None,
@@ -448,8 +450,9 @@ def __init__(
  diag_out = _spat_diag_out(self, w, 'yend', ml=True)
  else:
  self.title = ("MAXIMUM LIKELIHOOD SPATIAL LAG - REGIMES"+ " (METHOD = "+ method+ ")")
+
  if spat_impacts and slx_lags == 0:
- impacts = _summary_impacts(self, _spmultiplier(w, self.rho))
+ impacts = _summary_impacts(self, _spmultiplier(w, self.rho, method=spat_impacts), spat_impacts, regimes=True)
  try:
  diag_out += impacts
  except TypeError:
@@ -605,7 +608,7 @@ def ML_Lag_Regimes_Multi(
  if spat_diag and slx_lags == 1:
  results[r].other_mid += _spat_diag_out(results[r], None, 'yend', ml=True)
  if spat_impacts and slx_lags == 0:
- results[r].other_mid += _summary_impacts(results[r], _spmultiplier(results[r].w, results[r].rho))
+ results[r].other_mid += _summary_impacts(results[r], _spmultiplier(results[r].w, results[r].rho, method=spat_impacts), spat_impacts)
  counter += 1
  self.multi = results
  self.chow = REGI.Chow(self)
@@ -644,7 +647,7 @@ def _work(
  model.schwarz = DIAG.schwarz(reg=model)
  model.slx_lags = slx_lags
  model.w = w_r
- model.rho = model.betas[-1] 
+ model.rho = model.betas[-1]
  return model
 
 

spreg/output.py

@@ -5,10 +5,10 @@
 import textwrap as TW
 import numpy as np
 import pandas as pd
-import math
 from . import diagnostics as diagnostics
 from . import diagnostics_tsls as diagnostics_tsls
 from . import diagnostics_sp as diagnostics_sp
+from .sputils import _sp_effects
 
 __all__ = []
 
@@ -568,7 +568,7 @@ def _summary_iteration(reg):
 
  return txt
 
-def _summary_impacts(reg, spmult):
+def _summary_impacts(reg, spmult, spat_impacts, slx_lags=0, regimes=False):
  """
  Spatial direct, indirect and total effects in spatial lag model.
  Uses multipliers computed by sputils._spmultipliers.
@@ -577,21 +577,24 @@ def _summary_impacts(reg, spmult):
  ----------
  reg: spreg regression object
  spmult: spatial multipliers as a dictionary
+ spat_impacts: spatial impacts method as string
+ slx_lags: int, number of spatial lags of X in the model
+ regimes: boolean, True if regimes model
 
  Returns
  -------
  strings with direct, indirect and total effects
 
  """
- variables = reg.output[reg.output['var_type'].isin(['x', 'wx', 'yend']) & (reg.output.index != 0)]
+ variables = reg.output.query("var_type in ['x', 'yend'] and index != 0")
+ if regimes:
+ variables = variables[~variables['var_names'].str.endswith('_CONSTANT')]
  variables_index = variables.index
- m1 = spmult['ati']
- btot = m1 * reg.betas[variables_index]
- m2 = spmult['adi']
- bdir = m2 * reg.betas[variables_index]
- m3 = spmult['aii']
- bind = m3 * reg.betas[variables_index]
+
+ btot, bdir, bind = _sp_effects(reg, variables, spmult, slx_lags)
+
  strSummary = "\nSPATIAL LAG MODEL IMPACTS\n"
+ strSummary += "Impacts computed using the '" + spat_impacts + "' method.\n"
  strSummary += " Variable Direct Indirect Total\n"
  for i in range(len(variables)):
  strSummary += "%20s %12.4f %12.4f %12.4f\n" % (

spreg/sputils.py

@@ -307,8 +307,7 @@ def _spmultiplier(w, rho, method="simple", mtol=0.00000001):
  adii0 = np.sum(np.diag(invirw0))
  multipliers["adi"] = adii0 / n
  elif method == "power":
- wf = w.full()[0]
- ws3 = SP.csr_array(wf)
+ ws3 = w.to_sparse(fmt='csr') 
  rhop = rho
  ww = ws3
  pow = 1
@@ -329,6 +328,45 @@ def _spmultiplier(w, rho, method="simple", mtol=0.00000001):
  multipliers["aii"] = multipliers["ati"] - multipliers["adi"]
  return (multipliers)
 
+def _sp_effects(reg, variables, spmult, slx_lags=0):
+ """
+ Calculate spatial lag, direct and indirect effects
+ 
+ Attributes
+ ----------
+ reg : regression object
+ variables : chunk of self.output with variables to calculate effects
+ spmult : dictionary with spatial multipliers
+ slx_lags : number of SLX lags
+
+ Returns
+ -------
+ btot : total effects
+ bdir : direct effects
+ bind : indirect effects
+ """
+ variables_index = variables.index
+ m1 = spmult['ati']
+ btot = m1 * reg.betas[variables_index]
+ m2 = spmult['adi']
+ bdir = m2 * reg.betas[variables_index]
+
+ # Assumes all SLX effects are indirect effects. Needs revision by LA. 
+ if slx_lags > 0:
+ variables_wx = reg.output.query("var_type == 'wx'")
+ variables_wx_index = variables_wx.index
+
+ chunk_size = len(variables)
+ for i in range(slx_lags):
+ start_idx = i * chunk_size
+ end_idx = start_idx + chunk_size
+ chunk_indices = variables_wx_index[start_idx:end_idx]
+ btot += m1 * reg.betas[chunk_indices]
+ bind = btot - bdir
+ else:
+ m3 = spmult['aii']
+ bind = m3 * reg.betas[variables_index]
+ return btot, bdir, bind
 
 def _test():
  import doctest

spreg/twosls_sp.py

@@ -238,9 +238,11 @@ class GM_Lag(BaseGM_Lag):
  If True, then use n-k to estimate sigma^2. If False, use n.
  spat_diag : boolean
  If True, then compute Anselin-Kelejian test and Common Factor Hypothesis test (if applicable)
- spat_impacts : boolean
- If True, include average direct impact (ADI), average indirect impact (AII),
- and average total impact (ATI) in summary results
+ spat_impacts : string
+ Include average direct impact (ADI), average indirect impact (AII),
+ and average total impact (ATI) in summary results.
+ Options are 'simple', 'full', 'power', or None.
+ See sputils.spmultiplier for more information.
  vm : boolean
  If True, include variance-covariance matrix in summary
  results
@@ -505,7 +507,7 @@ def __init__(
  gwk=None,
  sig2n_k=False,
  spat_diag=True,
- spat_impacts=True,
+ spat_impacts="simple",
  vm=False,
  name_y=None,
  name_x=None,
@@ -585,10 +587,11 @@ def __init__(
  self.output['regime'], self.output['equation'] = (0, 0)
  self.other_top = _spat_pseudo_r2(self)
  diag_out = None
+
  if spat_diag:
  diag_out = _spat_diag_out(self, w, 'yend')
  if spat_impacts and slx_lags == 0:
- impacts = _summary_impacts(self, _spmultiplier(w, self.rho))
+ impacts = _summary_impacts(self, _spmultiplier(w, self.rho, method=spat_impacts), spat_impacts, slx_lags)
  try:
  diag_out += impacts
  except TypeError:

spreg/twosls_sp_regimes.py

@@ -88,9 +88,11 @@ class GM_Lag_Regimes(TSLS_Regimes, REGI.Regimes_Frame):
  matrix must have ones along the main diagonal.
  sig2n_k : boolean
  If True, then use n-k to estimate sigma^2. If False, use n.
- spat_impacts : boolean
- If True, include average direct impact (ADI), average indirect impact (AII),
- and average total impact (ATI) in summary results
+ spat_impacts : string
+ Include average direct impact (ADI), average indirect impact (AII),
+ and average total impact (ATI) in summary results.
+ Options are 'simple', 'full', 'power', or None.
+ See sputils.spmultiplier for more information.
  spat_diag : boolean
  If True, then compute Anselin-Kelejian test and Common Factor Hypothesis test (if applicable)
  vm : boolean
@@ -455,7 +457,7 @@ def __init__(
  gwk=None,
  sig2n_k=False,
  spat_diag=True,
- spat_impacts=True,
+ spat_impacts="simple",
  constant_regi="many",
  cols2regi="all",
  regime_lag_sep=False,
@@ -643,7 +645,7 @@ def __init__(
  if spat_diag:
  diag_out = _spat_diag_out(self, w, 'yend')
  if spat_impacts and slx_lags == 0:
- impacts = _summary_impacts(self, _spmultiplier(w, self.rho))
+ impacts = _summary_impacts(self, _spmultiplier(w, self.rho, method=spat_impacts), spat_impacts, regimes=True)
  try:
  diag_out += impacts
  except TypeError:
@@ -847,7 +849,7 @@ def GM_Lag_Regimes_Multi(
  if spat_diag:
  results[r].other_mid += _spat_diag_out(results[r], results[r].w, 'yend')
  if spat_impacts and slx_lags == 0:
- results[r].other_mid += _summary_impacts(results[r], _spmultiplier(results[r].w, results[r].rho))
+ results[r].other_mid += _summary_impacts(results[r], _spmultiplier(results[r].w, results[r].rho, method=spat_impacts), spat_impacts)
  counter += 1
  self.multi = results
  if robust == "hac":