[R-sig-Geo] SLX model for splm package in R

Tobias Rüttenauer ruettenauer at sowi.uni-kl.de
Fri Dec 22 11:18:42 CET 2017


Dear all,

though this is an old issue by now (but still one of the first google
results on the topic), here is some code answering one of the questions:
should you demean the spatial lag or lag the demeaned variable.

As far as I understand the following results, it does not matter. Both ways
produce identical results, which both conform to the LSDV approach.

Please find the example code below:


library(spdep)
library(splm)
library(texreg)

### Demean function
dm<-function(x, id){
  res<-x-ave(x,id,FUN=function(x) mean(x, na.rm=T))
  res
}

data(Produc, package="plm")
data(usaww)

usa.lw<-mat2listw(usaww)

Produc.pd<-pdata.frame(Produc, index=c("state", "year"))

Produc.pd$Wpcap<-slag(Produc.pd$pcap, usa.lw)
Produc.pd$Wpc<-slag(Produc.pd$pc, usa.lw)
Produc.pd$Wunemp<-slag(Produc.pd$unemp, usa.lw)

#### LSDV ####

lsdv.mod<-lm(gsp~pcap+pc+unemp + Wpcap+Wpc+Wunemp
             +state,
             data=Produc.pd)
summary(lsdv.mod)

#### Demean after lag ####

dal.mod<-plm(gsp~pcap+pc+unemp + Wpcap+Wpc+Wunemp,
             data=Produc.pd, effect="individual", model="within")
summary(dal.mod)

# Manual demeaning
Produc.pd[,4:ncol(Produc.pd)]<-apply(Produc.pd[,4:ncol(Produc.pd)], 2,
function(x) dm(x, Produc.pd$state))

dal2.mod<-lm(gsp~pcap+pc+unemp + Wpcap+Wpc+Wunemp,
             data=Produc.pd)
summary(dal2.mod)

#### Lag after demean ####

# Replace slags by demeaned slags
Produc.pd$Wpcap<-slag(Produc.pd$pcap, usa.lw)
Produc.pd$Wpc<-slag(Produc.pd$pc, usa.lw)
Produc.pd$Wunemp<-slag(Produc.pd$unemp, usa.lw)

lad.mod<-lm(gsp~pcap+pc+unemp + Wpcap+Wpc+Wunemp,
             data=Produc.pd)
summary(lad.mod)


screenreg(l=list(lsdv.mod, dal.mod, dal2.mod, lad.mod),
          omit.coef = "state", digits=6)



#### Test if both ways produce identical vectors ####

Produc.pd<-pdata.frame(Produc, index=c("state", "year"))

x<-Produc.pd$gsp

# lag after demean
dx<-dm(x, Produc.pd$state)
wdx<-slag(dx, usa.lw)

# deaman after lag
wx<-slag(x, usa.lw)
dwx<-dm(wx, Produc.pd$state)

# Test for equality
all.equal(wdx, dwx)


Best,
Tobi


Tobias Rüttenauer
TU Kaiserslautern
Erwin-Schrödinger-Straße 57
67663 Kaiserslautern

ruettenauer at sowi.uni-kl.de
Tel.: 0631 205 5785




> -----Ursprüngliche Nachricht-----
> Von: R-sig-Geo [mailto:r-sig-geo-bounces at r-project.org] Im Auftrag von
> Roger Bivand
> Gesendet: 30 August 2016 08:44
> An: dfamaral at usp.br
> Cc: r-sig-geo at r-project.org
> Betreff: Re: [R-sig-Geo] SLX model for splm package in R
> 
> On Tue, 30 Aug 2016, dfamaral at usp.br wrote:
> 
> > Dear all,
> >
> > Impacts in SLX model are calculated directly by parameters estimation:
> >
> > http://onlinelibrary.wiley.com/doi/10.1111/jors.12188/abstract
> >
> > Beta1 gives direct impact and Beta 2 the indirect:
> >
> > y = beta1*X + beta2*WX + e
> >
> 
> Yes, but you need the total impact beta1+beta2, but more work to infer
from
> that sum, which is what you need to assess the impact of X.
> 
> > I understand W must multiply T times the X matrix (T is the lenght of
> > time series for a panel data). And variables are demeaned after WX is
> > calculated.
> 
> Kronecker product of sparse matrices. Please contribute a patch for
> create_WX if desirable.
> 
> Please motivate the demean-after-lag concusion. Will this not affect the
> direct link between X and WX?
> 
> Roger
> 
> >
> > Daniel.
> >
> > ----- Mensagem original -----
> >
> >
> > De: "Roger Bivand" <Roger.Bivand at nhh.no>
> > Para: "Burcu Ozuduru" <bozuduru at gmail.com>
> > Cc: dfamaral at usp.br, r-sig-geo at r-project.org
> > Enviadas: Segunda-feira, 29 de Agosto de 2016 4:13:26
> > Assunto: Re: [R-sig-Geo] SLX model for splm package in R
> >
> > This is far too ad-hoc a discussion. For Spatial panel SLX, please
> > read
> > carefully:
> >
> > http://dx.doi.org/10.1016/j.spasta.2014.02.002
> > DOI 10.1007/s10109-015-0217-3
> > DOI: 10.1111/jors.12188
> >
> > for a view of the state of play. One question - should you demean the
> > X variables before taking spatial lags? This does not apply to
> > cross-section models for which spdep::create_WX() was written, for
> > obvious reasons. A second question - how do you handle the calculation
> > of total impacts in this setting? In cross-section models, this is
> > reasonably easy, by linear combination. But in spatial panel SLX
> > models, is it easy? What does the literature suggest?
> >
> > On Mon, 29 Aug 2016, Burcu Ozuduru wrote:
> >
> >> Hi,
> >>
> >> This is a great question because I faced the same problem a while back.
> >> What I did is to obtain the spatial weight matrix as (have spdep
> >> installed).
> >>
> >> S47matrix<-W1%*%S47
> >>
> >> and then put it in my data frame and equation. Btw, above W1 is a
> >> spatial weights matrix I obtained from a neighbor list
> >> (W1=as.matrix(as_dgRMatrix_listw(ccListw))) and S47 is the dependent
> >> variable.
> >>
> >> m2<-
> data.frame(S41,S46,S47,S65,S68,S47matrix,MAD10000,NQPDA10000,TPBt
> >> A10000,DivA10000,Lnk10000)
> >> summary(m2<-zeroinfl(S47~MAD10000+NQPDA10000+S47matrix,
> >> dist="negbin"))
> >>
> >> I also wanted to ask whether I was on the right track. Professor
> >> Bivand, if this is a common problem could you direct us to a relevant
> >> page -if available and if at all possible- please? Your advice would
> >> be sincerely appreciated.
> >
> > This is an inappropriate question for this thread (SLX spatial panel
> > models). You must refer to the literature, and understand what you are
> > doing:
> >
> > zeroinfl(S47~MAD10000+NQPDA10000+S47matrix, dist="negbin")
> >
> > is completely wrong, as S47matrix is Wy, the spatial lag of the
> > dependent variable, and this spatial lag model of response S47 cannot
> > be estimated using pscl::zeroinfl (do say which package provides the
> > functions you mention). Your only recourse is to write your own
> > Bayesian model. It is even extremely likely that this model doesn't
> > make any substantive sense, unless you have a good understanding about
> > why S47 at i (for all i) should be simultaneously caused by the full
> > set of values of S47. A model in the error might make sense (a zero
> > inflated negbin model with an additive mrf random effect (maybe
> R2BayesX 10.18637/jss.v014.i11).
> >
> > An SLX zero inflated negbin model could maybe apply, in which case
> > pscl::zeroinfl might be applicable, but then you'd be including the
> > lags of MAD10000 and NQPDA10000, not the response, and inference on
> > total impacts of these variables could be challenging rather than
> > obvious (MCMC summaries of the summed samples of the coefficients
> > would work, but here no MCMC I think - R2BayesX might help).
> >
> > Most of these issues are far more complex than students or their
> > supervisors understand, do ask an experienced applied statistician
> > (I'm not an applied statistician, don't rely on my suggestions, you
> > yourself are responsible for your results).
> >
> > Hope this clarifies,
> >
> > Roger
> >
> >>
> >> Kind Regards-
> >>
> >> Thanks-
> >>
> >> Burcu
> >>
> >> ---
> >>
> >> Burcu H. Ozuduru,
> >> Gazi University, Faculty of Architecture, Department of City and
> >> Regional Planning Maltepe-Ankara Turkey
> >> t: +90-312-5823701
> >> e: bozuduru at gazi.edu.tr
> >> http://websitem.gazi.edu.tr/bozuduru
> >>
> >>
> >> On Mon, Aug 29, 2016 at 12:42 AM, Roger Bivand
> <Roger.Bivand at nhh.no> wrote:
> >>
> >>> On Sun, 28 Aug 2016, dfamaral at usp.br wrote:
> >>>
> >>> Dear colleagues,
> >>>>
> >>>> Does anyone know how to implement Spatial Lag of X (SLX) model for
> >>>> panel data using splm package? I tried the following to create a WX
> >>>> matrix and then apply Fixed and Random Effects:
> >>>>
> >>>> tb <- read.csv2("panel.csv", header = TRUE, sep = ";", dec = ".")
> >>>> time <- length(unique(tb$year)) X <- as.matrix(tb[ ,3:20])
> >>>>
> >>>
> >>> class(X) is not "pseries", is it?
> >>>
> >>> WX <- slag(X, listw, maxlag = 1)
> >>>>
> >>>> I see that slag is applicable for vectors of pseries, but this
command
> >>>> led to the error message: "Error in UseMethod("slag"): no applicable
> method
> >>>> for 'slag' applied to an object of class "c('matrix', 'double',
'numeric')"
> >>>>
> >>>
> >>> slag.default() takes a vector, slag.pseries() a "pseries" object.
> >>>
> >>>
> >>>> Is there any way to implement this SLX model for panel series?
> >>>>
> >>>
> >>> Did you look at spdep::create_WX?
> >>>
> >>>
> >>>> Regards,
> >>>>
> >>>> Daniel.
> >>>>
> >>>> [[alternative HTML version deleted]]
> >>>>
> >>>>
> >>> Please do not post HTML - like insects, HTML carries infectious
payloads
> >>> and wastes bandwidth and server processing capacity.
> >>>
> >>> _______________________________________________
> >>>> R-sig-Geo mailing list
> >>>> R-sig-Geo at r-project.org
> >>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >>>>
> >>>>
> >>> --
> >>> Roger Bivand
> >>> Department of Economics, Norwegian School of Economics,
> >>> Helleveien 30, N-5045 Bergen, Norway.
> >>> voice: +47 55 95 93 55; fax +47 55 95 91 00
> >>> e-mail: Roger.Bivand at nhh.no
> >>> http://orcid.org/0000-0003-2392-6140
> >>> https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en
> >>> http://depsy.org/person/434412
> >>>
> >>> _______________________________________________
> >>> R-sig-Geo mailing list
> >>> R-sig-Geo at r-project.org
> >>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >>>
> >>
> >
> >
> 
> --
> Roger Bivand
> Department of Economics, Norwegian School of Economics,
> Helleveien 30, N-5045 Bergen, Norway.
> voice: +47 55 95 93 55; fax +47 55 95 91 00
> e-mail: Roger.Bivand at nhh.no
> http://orcid.org/0000-0003-2392-6140
> https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en
> http://depsy.org/person/434412
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo



More information about the R-sig-Geo mailing list