[R-sig-Geo] Spatial lags using gridded data

Roger Bivand Roger.Bivand at nhh.no
Thu Aug 21 12:15:28 CEST 2014


On Thu, 21 Aug 2014, Steve Pickering wrote:

> Hello, all,
>
> Sorry for the silly question.
>
> I'm struggling to apply spatial lags to a data grid.  I've been
> working through the notes on Luc Anselin's spdep, but just can't
> figure out how to apply it to gridded data.
>
> Consider the following very simple 5 x 5 grid:
>
> ids <- c(1:25)
> lats <- c(30,30,30,30,30,31,31,31,31,31,32,32,32,32,32,33,33,33,33,33,34,34,34,34,34)
> lngs <- c(10,11,12,13,14,10,11,12,13,14,10,11,12,13,14,10,11,12,13,14,10,11,12,13,14)
> conflict <- c(0,0,0,0,0,0,0,0,1,0,0,0,1,1,0,0,1,0,0,0,0,0,0,0,0)
> conflictLag <- c(0,0,0.125,0.125,0.125,0,0.125,0.375,0.25,0.25,0.125,0.25,0.375,0.25,0.25,0.125,0.125,0.375,0.25,0.125,0.125,0.125,0.125,0,0)
> pop <- sample(100, 25)
>
> newGrid = data.frame(ids, lats, lngs, conflict, conflictLag, pop)
>
> Here, I have a conflict dummy, to determine whether there was some
> sort of conflict or war in that grid cell.  I also have a population
> variable (for the purposes of this exercise, determined at random).
>
> If I wanted to run a regression on this without including spatial
> lags, I could do something very simple, like this:
>
> model1 = lm(newGrid$conflict ~ newGrid$pop)
> summary(model1)

Thanks for the toy example:

set.seed(1)
pop <- sample(100, 25)
newGrid = data.frame(ids, lats, lngs, conflict, conflictLag, pop)
model1 = lm(conflict ~ pop, data=newGrid)
summary(model1)

Then make a spatial object:

library(sp)
coordinates(newGrid) <- c("lngs", "lats")
gridded(newGrid) <- TRUE
image(newGrid, "conflict")

Create a list of neighbours (here ignoring the possibility that the grid 
is in geographical, not projected coordinates:

library(spdep)
nb <- dnearneigh(newGrid, 0, sqrt(2))
summary(nb)
plot(nb, coordinates(newGrid), add=TRUE)

Choose a model and fitting function:

model2 <- spautolm(conflict ~ pop, data=newGrid, listw=nb2listw(nb,
  style="W"), family="SAR")
summary(model2)

but note that conflict is not a continuous variable, so:

model1a = glm(conflict ~ pop, data=newGrid, family=binomial(link=probit))
summary(model1a)

and

library(spatialprobit)
W <- as(as_dgRMatrix_listw(nb2listw(nb, style="W")), "CsparseMatrix")
model2a <- semprobit(conflict ~ pop, W=W, data=newGrid)
summary(model2a)

which uses MCMC, so is slow. Alternatively, but for the spatially 
lagged response, not for the spatially lagged error term:

model2b <- sarprobit(conflict ~ pop, W=W, data=newGrid)
summary(model2b)

This object includes effects/impacts (see LeSage & Pace 2009, or Ward & 
Gleditsch 2008).

library(McSpatial)
model2c <- spprobit(conflict ~ pop, wmat=nb2mat(nb), data=newGrid)
model2d <- spprobitml(conflict ~ pop, wmat=nb2mat(nb), data=newGrid)

use GMM and ML respectively and do not provide impacts. A recent paper in 
Journal of Regional Science compares spatial probit estimation approaches. 
In addition, the tally of 0 conflicts is 4 of 25, so your data may well 
also be zero-inflated.

Spatial misspecification is possibly only a minor problem, really.

Hope this helps,

Roger


>
> However, I know that I am supposed to include a spatial lag.  As you
> can see in the code, I have also created a spatial lag variable based
> on whether there was conflict in adjacent cells (Queen's move, hence a
> fraction of eight).
>
> What I'm struggling with is the correct way to include the spatial lag
> in my model.
>
> Thanks, and sorry for the simple question.
>
> - Steve.
>
> _______________________________________________
> 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



More information about the R-sig-Geo mailing list