**for**

__archive__**R-sig-geo**, not a forum. It is not possible to post through Nabble - you may not start a new thread nor follow up an existing thread. If you wish to post, but are not subscribed to the list through the list homepage, subscribe first

**(through the list homepage, not via Nabble)**and post from your subscribed email address. Until 2015-06-20, subscribers could post through Nabble, but policy has been changed as too many non-subscribers misunderstood the interface.

*Updated:*1 hour 42 min ago

### Correlation between a nominal and a continious rasters

Hi all

I need to find how much two *raster* variables are correlated. One of these

variables is categorical (*nominal* such as land-cover) and the other one

is *continuous*.

I also need to do the same analysis for a *Binary raster layer* (0,1)

and a *continuous

raster layer*. Could you please give me some bits of advice about the best

methods in this regard in R?

Excuse me if you find my question so elementary because I’m new in R.

Best regards

Iman

[[alternative HTML version deleted]]

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### layer to change nodata values

> I downloaded global water surface layer from Pekel et al., 2016. The data (band values) has 0 to 255. 0 = 0 water, 1 = 1%, 2= 2%….100=100% water and after then 255 (which is no data according to paper). I want to keep this values from 0 to 100 only in the raster and want to change values 255 into nodata. How can I make r code? I have following code but I am not sure it is good not?

a <- raster(“waterlayer.tif)

max <- 100

a[a >= max] <- NA

Hari Sharma

Taiwan

[[alternative HTML version deleted]]

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: How to test spatial dependence in errorsarlm

> So testing for spatial dependence on the residuals by means of the

> lm.LMtests (option LMerr, the only that works with residuals) is wrong,

> isn't it? I had read in some forum that this was a posible way to test it...

There was some speculation that these tests might be used on an lm fit of

the (I - \lambda W) y ~ (I - \lambda W) X model. These speculations have

never been checked rigorously, so nobody knows whether they are of any

value, probably not, and certainly we know nothing of their inferential

basis.

>

> In my case the Moran test and the LM tests (both LMerr and LMlag, and

> also their robust versions) are strongly rejected (p-values between

> 4.307e-06 and 2.2e-16). As the rejection is stronger for the spatial

> error model, my suspicion was that this could be the best model to

> capture the spatial dependence (in fact the log-likelihood is bigger for

> the spatial error model, and the AIC lower). However, how can I know

> whether the spatial error model is a good option if I cannot test the

> absence of spatial dependence in the residuals?

You by definition fix X and W as known, exogeneous, quantities. If X (and

functional forms) and/or W do not meet these requirements, there is little

good guidance. It depends on what you want: predict house prices where

they are not observed; estimate \beta values; estimate the impact of a

unit change in an X variable on house prices (y); whatever. A best-fit

model suggests that you want to predict, but isn't necessary for impacts

or betas if you trust X (and its functional forms) and W.

> And how can I know, as you suspect, whether I have a

> misspecification problem?

See a good deal of work by Daniel McMillen on these issues.

> Moreover, I also estimated the Durbin model, and in this case the LM

> test on the residuals suggests no spatial dependence (for the spatial

> lag model I get the opposite conclusion), but due to the nature of my

> regression I don't think that this model is suitable (the regressors are

> characteristics of houses such as size, number of rooms, etc).

This suggests that SLX or SDEM (see LeSage 2014 and SLX articles for a

discussion) may address many of the issues of spatial autocorrelation by

including (selected) WX. The Durbin versions of spdep functions do not

(yet) let you choose which WX to include, always including all X - this is

on my medium-term to-do list.

Roger

>

> Thanks a lot for your time.

> Best

> Javi

>

> -----Mensaje original-----

> De: Roger Bivand [mailto:[hidden email]]

> Enviado el: domingo, 13 de agosto de 2017 12:45

> Para: Javier García

> CC: [hidden email]

> Asunto: Re: [R-sig-Geo] How to test spatial dependence in errorsarlm

>

> On Sun, 13 Aug 2017, Javier García wrote:

>

>> Hello everybody:

>>

>>

>>

>> I have estimated a spatial error model and now I would like to test

>> whether that model has really ?deleted? the spatial dependence. For

>> the spatial lag model and for the Durbin model the function lagsarlm

>> gives the LM test for residual autocorrelation test value, but the

> function errorsarlm does not.

>> Does anyone know how to do it in R?

>>

>

> As you should be aware from the literature, the only LM test that has been

> written (the maths) is a test for residual error autocorrelation for spatial

> lag models. Doing it in R will not help until someone (you?) does the maths.

> Computing a value is easy, but knowing what to infer from it is hard. By

> definition, if your model is well-specified, the residual autocorrelation is

> fully captured by its coefficient. I suspect that your model suffers from

> mis-specification problems.

>

> Roger

>

>>

>>

>> Thanks a lot in advance.

>>

>>

>>

>> Javi

>>

>>

>>

>>

>>

>>

>> JAVIER GARCÍA

>>

>>

>>

>> Departamento de Economía Aplicada III (Econometría y Estadística)

>>

>> Facultad de Economía y Empresa (Sección Sarriko)

>> Avda. Lehendakari Aguirre 83

>>

>> 48015 BILBAO

>> T.: +34 601 7126 F.: +34 601 3754

>> <http://www.ehu.es/> www.ehu.es

>>

>>

> http://www.unibertsitate-hedakuntza.ehu.es/p268-content/es/contenidos/inform

>>

> acion/manual_id_corp/es_manual/images/firma_email_upv_euskampus_bilingue.gif

>>

>>

>>

>>

>>

>>

>

> --

Roger Bivand

Department of Economics, Norwegian School of Economics,

Helleveien 30, N-5045 Bergen, Norway.

voice: +47 55 95 93 55; e-mail: [hidden email]

Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html

http://orcid.org/0000-0003-2392-6140

https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo Roger Bivand

Department of Economics

Norwegian School of Economics

Helleveien 30

N-5045 Bergen, Norway

### Re: How to test spatial dependence in errorsarlm

lm.LMtests (option LMerr, the only that works with residuals) is wrong,

isn't it? I had read in some forum that this was a posible way to test it...

In my case the Moran test and the LM tests (both LMerr and LMlag, and also

their robust versions) are strongly rejected (p-values between 4.307e-06

and 2.2e-16). As the rejection is stronger for the spatial error model, my

suspicion was that this could be the best model to capture the spatial

dependence (in fact the log-likelihood is bigger for the spatial error

model, and the AIC lower). However, how can I know whether the spatial error

model is a good option if I cannot test the absence of spatial dependence in

the residuals? And how can I know, as you suspect, whether I have a

misspecification problem? Moreover, I also estimated the Durbin model, and

in this case the LM test on the residuals suggests no spatial dependence

(for the spatial lag model I get the opposite conclusion), but due to the

nature of my regression I don't think that this model is suitable (the

regressors are characteristics of houses such as size, number of rooms,

etc).

Thanks a lot for your time.

Best

Javi

-----Mensaje original-----

De: Roger Bivand [mailto:[hidden email]]

Enviado el: domingo, 13 de agosto de 2017 12:45

Para: Javier García

CC: [hidden email]

Asunto: Re: [R-sig-Geo] How to test spatial dependence in errorsarlm

On Sun, 13 Aug 2017, Javier García wrote:

> Hello everybody:

>

>

>

> I have estimated a spatial error model and now I would like to test

> whether that model has really ?deleted? the spatial dependence. For

> the spatial lag model and for the Durbin model the function lagsarlm

> gives the LM test for residual autocorrelation test value, but the

function errorsarlm does not.

> Does anyone know how to do it in R?

>

As you should be aware from the literature, the only LM test that has been

written (the maths) is a test for residual error autocorrelation for spatial

lag models. Doing it in R will not help until someone (you?) does the maths.

Computing a value is easy, but knowing what to infer from it is hard. By

definition, if your model is well-specified, the residual autocorrelation is

fully captured by its coefficient. I suspect that your model suffers from

mis-specification problems.

Roger

>

>

> Thanks a lot in advance.

>

>

>

> Javi

>

>

>

>

>

>

> JAVIER GARCÍA

>

>

>

> Departamento de Economía Aplicada III (Econometría y Estadística)

>

> Facultad de Economía y Empresa (Sección Sarriko)

> Avda. Lehendakari Aguirre 83

>

> 48015 BILBAO

> T.: +34 601 7126 F.: +34 601 3754

> <http://www.ehu.es/> www.ehu.es

>

> http://www.unibertsitate-hedakuntza.ehu.es/p268-content/es/contenidos/inform

>

acion/manual_id_corp/es_manual/images/firma_email_upv_euskampus_bilingue.gif

>

>

>

>

>

>

--

Roger Bivand

Department of Economics, Norwegian School of Economics,

Helleveien 30, N-5045 Bergen, Norway.

voice: +47 55 95 93 55; e-mail: [hidden email]

Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html

http://orcid.org/0000-0003-2392-6140

https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: Heteroscedasticity in Spatial Error Model

> Hello again:

>

>

>

> Please, could anyone tell me how to estimate robust standard errors for a

> spatial error model? The residuals of my model show heteroscedasticity

> evidence, but the functions I have looked at only work with lm type objects.

>

The documented approaches use GM rather than ML for fitting - see the

sphet package and https://www.jstatsoft.org/index.php/jss/issue/view/v063.

You should really try to remove the sources of model mis-specification

instead of spreading coefficient standard errors by guesswork. The may

stem from MAUP, missing covariates and/or wrong functional forms.

Roger

>

>

> Thanks a lot in advance.

>

>

>

> Javi

>

>

>

>

>

>

> JAVIER GARCÍA

>

>

>

> Departamento de Economía Aplicada III (Econometría y Estadística)

>

> Facultad de Economía y Empresa (Sección Sarriko)

> Avda. Lehendakari Aguirre 83

>

> 48015 BILBAO

> T.: +34 601 7126 F.: +34 601 3754

> <http://www.ehu.es/> www.ehu.es

>

> http://www.unibertsitate-hedakuntza.ehu.es/p268-content/es/contenidos/inform

> acion/manual_id_corp/es_manual/images/firma_email_upv_euskampus_bilingue.gif

>

>

>

>

>

> --

Roger Bivand

Department of Economics, Norwegian School of Economics,

Helleveien 30, N-5045 Bergen, Norway.

voice: +47 55 95 93 55; e-mail: [hidden email]

Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html

http://orcid.org/0000-0003-2392-6140

https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo Roger Bivand

Department of Economics

Norwegian School of Economics

Helleveien 30

N-5045 Bergen, Norway

### Re: How to test spatial dependence in errorsarlm

> Hello everybody:

>

>

>

> I have estimated a spatial error model and now I would like to test whether

> that model has really ?deleted? the spatial dependence. For the spatial lag

> model and for the Durbin model the function lagsarlm gives the LM test for

> residual autocorrelation test value, but the function errorsarlm does not.

> Does anyone know how to do it in R?

> As you should be aware from the literature, the only LM test that has been

written (the maths) is a test for residual error autocorrelation for

spatial lag models. Doing it in R will not help until someone (you?) does

the maths. Computing a value is easy, but knowing what to infer from it is

hard. By definition, if your model is well-specified, the residual

autocorrelation is fully captured by its coefficient. I suspect that your

model suffers from mis-specification problems.

Roger

>

>

> Thanks a lot in advance.

>

>

>

> Javi

>

>

>

>

>

>

> JAVIER GARCÍA

>

>

>

> Departamento de Economía Aplicada III (Econometría y Estadística)

>

> Facultad de Economía y Empresa (Sección Sarriko)

> Avda. Lehendakari Aguirre 83

>

> 48015 BILBAO

> T.: +34 601 7126 F.: +34 601 3754

> <http://www.ehu.es/> www.ehu.es

>

> http://www.unibertsitate-hedakuntza.ehu.es/p268-content/es/contenidos/inform

> acion/manual_id_corp/es_manual/images/firma_email_upv_euskampus_bilingue.gif

>

>

>

>

>

> --

Roger Bivand

Department of Economics, Norwegian School of Economics,

Helleveien 30, N-5045 Bergen, Norway.

voice: +47 55 95 93 55; e-mail: [hidden email]

Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html

http://orcid.org/0000-0003-2392-6140

https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo Roger Bivand

Department of Economics

Norwegian School of Economics

Helleveien 30

N-5045 Bergen, Norway

### Heteroscedasticity in Spatial Error Model

Hello again:

Please, could anyone tell me how to estimate robust standard errors for a spatial error model? The residuals of my model show heteroscedasticity evidence, but the functions I have looked at only work with lm type objects.

Thanks a lot in advance.

Javi

**JAVIER GARCÍA **

** **

**Departamento de Economía Aplicada III (Econometría y Estadística)**

**Facultad de Economía y Empresa (Sección Sarriko)**

Avda. Lehendakari Aguirre 83

48015 BILBAO**T.: +34 601 7126** **F.: +34 601 3754****www.ehu.es**

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### How to test spatial dependence in errorsarlm

Hello everybody:

I have estimated a spatial error model and now I would like to test whether that model has really “deleted” the spatial dependence. For the spatial lag model and for the Durbin model the function lagsarlm gives the LM test for residual autocorrelation test value, but the function errorsarlm does not. Does anyone know how to do it in R?

Thanks a lot in advance.

Javi

**JAVIER GARCÍA **

** **

**Departamento de Economía Aplicada III (Econometría y Estadística)**

**Facultad de Economía y Empresa (Sección Sarriko)**

Avda. Lehendakari Aguirre 83

48015 BILBAO**T.: +34 601 7126** **F.: +34 601 3754****www.ehu.es**

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Regression with rasters using calc::raster

I'm trying to fit the Firth's Penalized MLE GLM implemented in logistf

package to a set of rasters but I'm facing errors and problems I couldn't

realize until now.

# Lets generate 2 rasters reproducing what I'm facing

r <- raster(nrow=10, ncol=10)

# binary response raster-variable with 9 bands

s1 <- lapply(1:9, function(i) setValues(r, sample(0:1,ncell(r),replace =

T)))

s1 <- stack(s1)

# one explanatory raster-variable

val <- sample(0:60,ncell(r),replace = T)

s2 <- raster(nrow=10, ncol=10,vals=val)

plot(s1)

plot(s2)

# a second explanatory variable. Nine values

exp_2 <- c(27.00,30.02,31.07,32.72,33.73,35.12,35.65,36.06,38.32)

Now, I want to fit a model using Firth's Penalized MLE GLM implemented in

logistf (i have reasons for this not reproduced here) using calc from

raster package. That's where the mystery lives.

The rationale is each cell in:

s1/layer1 ~ 27.00 + corresponding cell in s2 + 27.00:corresponding cell in

s2

s1/layer2 ~ 30.02 + corresponding cell in s2 + 30.02:corresponding cell in

s2

s1/layer3 ~ 31.07 + corresponding cell in s2 + 31.07:corresponding cell in

s2

... and so on for all 9 bands of my response raster-variable, which are

paired with values from exp_2.

# I tried something like this:

fun <- function(x) { logistf(x ~ exp_2 + s2 + exp_2:s2)$coefficients }

coefs <- calc(s1,fun)

But it was clear it wouldn't work. The tricky part is to tell R I want each

value of exp_2 to be used with each rasterlayer of s1 for this model.

Any idea would be appreciated. Ideas?

Thanks in advance

*Jefferson Ferreira-Ferreira, **PhD (abd)*

*Geographer*

*Ecosystem Dynamics Observatory <http://tscanada.wix.com/ecodyn> -

EcoDyn/UNESP*

*Department of **Geography *

*Institute of Geosciences and Exact Sciences** (IGCE)*

*São Paulo State University (UNESP)*

*Rio Claro, SP - Brazil*

[[alternative HTML version deleted]]

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### $columns that I added to a data frame don't have attr(*, "label")

(that I made)

have attr() like the $columns that came with the data

$ CLUSTER :Classes 'labelled', 'integer' atomic [1:90435] 4 4 4 4

4 4 4 4 4 4 ...

.. ..- attr(*, "label")= Named chr "Socio-economic level of locality"

.. .. ..- attr(*, "names")= chr "CLUSTER"

$ size : chr "largest" "largest" "largest" "largest" ...

$ stock : num 596 596 596 596 596 ...

[[alternative HTML version deleted]]

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: help with plotKML::kml()

I think I set the color of the polygon line "white" which is the KML

default, so you can only change the color of polygons (fill); see:

https://github.com/cran/plotKML/blob/master/R/layer.SpatialPolygons.R#L112

To change colors per polygon you should run:

data(eberg_zones)

legend = get("colour_scale_factor", envir = plotKML.opts)

kml_open("eberg_zones.kml")

kml_layer(eberg_zones, colour=ZONES, colour_scale=legend)

kml_close("eberg_zones.kml")

See also: http://gsif.isric.org/doku.php/wiki:tutorial_plotkml

Polygon outlines

(https://developers.google.com/kml/documentation/kmlreference#polystyle)

could also be changed but not via plotKML.

HTH,

On 07-08-17 17:46, MacQueen, Don wrote:

> Hi Tomislav,

>

> I am using the kml() function in the plotKML package to write a SpatialPolygonsDataFrame object to a kml file, and would like to know if it is possible to specify the fill color and border color separately. I have been looking at documentation and been unable to find a way. I would appreciate any suggestions.

>

>

> Suppose "foo" is a SpatialPolygonsDataFrame object with a lat/long coordinate reference system. Then in R I can do, for example,

>

> plot(foo, col='yellow', border='red', lwd=2)

>

> I would like to write a KML file that will display the polygon in Google Earth in a similar manner.

>

> I have found that

> plotKML::kml(foo, colour='yellow')

> controls the fill color, but I can't find how to set the border color.

>

> Ideally, I'd also like to be able specify border only, no fill, and also line widths for the border.

>

> For example, the kml file might have in it something like:

> <Style id="blue">

> <PolyStyle>

> <color>4DFF0000</color>

> <fill>1</fill>

> <outline>1</outline>

> </PolyStyle>

> <LineStyle>

> <color>FFFF0000</color>

> <width>2</width>

> </LineStyle>

> </Style>

>

>

> Thanks

> -Don

>

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: How to draw the same legend (one legend) for the multiple spatial figures?

probably the most straightforward way to combine multiple trellis plots

into one is via latticeExtra::c.trellis(). Be aware, however, that the

desired range of z-values (from which your legend will ultimately be

created) needs to be assigned manually to each sub-plot using 'at'. Here

is a minimal example based on the meuse.grid 'SpatialPixelsDataFrame'

from sp, which works the same eg for 'Raster*' objects.

-----

## sample data

library(sp)

data("meuse.grid")

gridded(meuse.grid) = ~ x + y

meuse.grid$dist2 = meuse.grid$dist^2

p1 = spplot(meuse.grid, zcol = "dist", at = seq(0, 1, .01),

sp.layout = list("sp.text", c(179000, 333250), "a) dist"),

colorkey = list(height = .5), scales = list(draw = TRUE))

p2 = spplot(meuse.grid, zcol = "dist2", at = seq(0, 1, .01),

sp.layout = list("sp.text", c(179000, 333250), "b) dist2"))

## combine plots using latticeExtra::c.trellis

update(c(p1, p2), layout = c(1, 2), as.table = TRUE) # 1 column, 2 rows

update(c(p1, p2)) # 1 row, 2 columns

-----

In order to combine numerous plots into one (eg stored in a 'list'), you

can take inspiration from the Reduce()-based approach in

Orcs::latticeCombineGrid() (see

https://github.com/fdetsch/Orcs/blob/master/R/latticeCombineGrid.R;

credit goes to Tim Appelhans), among others.

Best,

Florian

On 06.08.2017 14:13, Bede-Fazekas Ákos wrote:

> Dear Brandon and Wei,

> I don't know the answer to your question but a somewhat similar

> solution is when you plot a separate legend using lattice::draw.key()

> and disable sp::spplot()'s built-in legend using argument "auto.key =

> FALSE".

> An example:

>

> # loading libraries

> library(lattice)

> library(sp)

> library(gridExtra)

> library(grid)

>

> # drawing legend

> no_of_categories <- 10

> cutpoints <- seq(from = 0 - 0.0001, to = 1 + 0.0001, length.out =

> no_of_categories + 1)

> legend_text <- apply(X = cbind(round(cutpoints [1:no_of_categories],

> 1), round(cutpoints [2:(no_of_categories + 1)], 1)), MARGIN = 1, FUN =

> function(number) {return(paste(format(x = number, digits = 2),

> collapse = " - "))})

> legend_colors <- colorRampPalette(c("red", "orange", "yellow",

> "lightgreen", "darkgreen"))(no_of_categories)

> legend <- draw.key(key = list(reverse.rows = TRUE, space = "right",

> rectangles = list(col = legend_colors, border = FALSE), text =

> list(legend_text)), draw = FALSE)

> dev.off()

>

> # drawing maps

> maps <- list()

> for (column_name in colnames(sp_object@data)) {

> map <- spplot(obj = sp_object, zcol = column_name, auto.key =

> FALSE, col.regions = legend_colors, cuts = cutpoints)

> maps <- append(maps , list(map))

> }

> maps <- append(maps , list(legend))

>

> # plotting the map-legend composite to a png file

> layout <- rbind(c(1,2), c(3,4), c(5,6)) # let's say we have 5 maps and

> a legend

> composite <- do.call(arrangeGrob, c(maps, list(layout_matrix =

> layout), list(widths = c(1, 1))))

> png(width = 1000, height = 1000, filename = "map.png"))

> grid.newpage()

> grid.draw(composite)

> dev.off()

>

> Hope this helps,

> Ákos Bede-Fazekas

> Hungarian Academy of Sciences

>

> 2017.08.06. 11:53 keltezéssel, Brandon Payne írta:

>> How to draw the same legend (one legend) for the

>> multiple spatial figures?

>>

>> I would put the legend next to the upper-most plot and

>> wrap the whole thing in a single

>>

>> /figure

>>

>> so that it would be more obvious that the same legend applied to all.

>>

>> _______________________________________________

>> R-sig-Geo mailing list

>> [hidden email]

>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

>>

>

> _______________________________________________

> R-sig-Geo mailing list

> [hidden email]

> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

>

--

Dr. Florian Detsch

Environmental Informatics

Department of Geography

Philipps-Universität Marburg

Deutschhausstraße 12

35032 (parcel post: 35037) Marburg, Germany

Phone: +49 (0) 6421 28-25323

Web: http://www.uni-marburg.de/fb19/fachgebiete/umweltinformatik/detschf/index.html

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: How to draw the same legend (one legend) for the multiple spatial figures?

I don't know the answer to your question but a somewhat similar solution

is when you plot a separate legend using lattice::draw.key() and disable

sp::spplot()'s built-in legend using argument "auto.key = FALSE".

An example:

# loading libraries

library(lattice)

library(sp)

library(gridExtra)

library(grid)

# drawing legend

no_of_categories <- 10

cutpoints <- seq(from = 0 - 0.0001, to = 1 + 0.0001, length.out =

no_of_categories + 1)

legend_text <- apply(X = cbind(round(cutpoints [1:no_of_categories], 1),

round(cutpoints [2:(no_of_categories + 1)], 1)), MARGIN = 1, FUN =

function(number) {return(paste(format(x = number, digits = 2), collapse

= " - "))})

legend_colors <- colorRampPalette(c("red", "orange", "yellow",

"lightgreen", "darkgreen"))(no_of_categories)

legend <- draw.key(key = list(reverse.rows = TRUE, space = "right",

rectangles = list(col = legend_colors, border = FALSE), text =

list(legend_text)), draw = FALSE)

dev.off()

# drawing maps

maps <- list()

for (column_name in colnames(sp_object@data)) {

map <- spplot(obj = sp_object, zcol = column_name, auto.key =

FALSE, col.regions = legend_colors, cuts = cutpoints)

maps <- append(maps , list(map))

}

maps <- append(maps , list(legend))

# plotting the map-legend composite to a png file

layout <- rbind(c(1,2), c(3,4), c(5,6)) # let's say we have 5 maps and a

legend

composite <- do.call(arrangeGrob, c(maps, list(layout_matrix = layout),

list(widths = c(1, 1))))

png(width = 1000, height = 1000, filename = "map.png"))

grid.newpage()

grid.draw(composite)

dev.off()

Hope this helps,

Ákos Bede-Fazekas

Hungarian Academy of Sciences

2017.08.06. 11:53 keltezéssel, Brandon Payne írta:

> How to draw the same legend (one legend) for the

> multiple spatial figures?

>

> I would put the legend next to the upper-most plot and

> wrap the whole thing in a single

>

> /figure

>

> so that it would be more obvious that the same legend applied to all.

>

> _______________________________________________

> R-sig-Geo mailing list

> [hidden email]

> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

>

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: bivariate spatial correlation in R

you're correct in pointing out that the code would much more work to be

incorporated in the spdep library. This is something that Rogerio and I are

considering to do in the next few months as both of us find more time in

our agendas.

Thank you again for the discussion. The amount of work and attention you

put into this community is exceptional for someone as busy as you are and

this is much appreciated!

best,

Rafa

On Thu, Aug 3, 2017 at 10:13 AM, Roger Bivand <[hidden email]> wrote:

> On Wed, 2 Aug 2017, Rafael Pereira wrote:

>

> Roger,

>>

>>

>> Thank you for your response.

>>

>>

>> I recognize the data is not ideal and the analysis has limitations because

>> of the lack of information on population displacements that might have

>> occurred over the years. Nonetheless, there are plenty of data +

>> literature

>> showing how the spatial distribution of income classes and land use

>> patterns is fairly stable over time in this city, particularly for a short

>> timescales like in this analysis. Having said this, I believe these two

>> questions (1) what socioeconomic classes have gained more accessibility?

>> and (2) “were wealthier areas in 2010 able to attract more changes to

>> accessibility?” in the end ask the same thing but with different

>> phrasings,

>> though your phrasing (2) is more precise/correct.

>>

>>

>>

>> On the more technical discussion, I see your point that spatial AND

>> temporal correlation in my data would make permutation bootstrap

>> inappropriate to generate significance levels, thus making bivariate

>> Moran’s I biased. Thank you very much for the clarifications! This has

>> been

>> very helpful and I will have a closer look at which spatial regression

>> models are more appropriate for my analysis.

>>

>>

>> On a side note, do you think the function to calculate bivariate Moran’s I

>> is correct? And could it be incorporated in the next version of spdep? If

>> so, please give credit to Rogério Barbosa, the researcher who proposed the

>> code in Stack Overflow.

>>

>

> Perhaps, but SO is usually ephemeral (nobody takes responsibility for

> documenting code and fixing bugs found later). I don't see any tests,

> documentation or accommodation of what spdep expects for edge cases - the

> function as it stands would need a lot of work to protect users from

> obvious blunders. There are no references to literature, nor to proven test

> cases which this implementation should match. We have an implementation of

> Lee (2001), but this is not the same, right? Which article gives the formal

> statistical development of the bivariate local Moran's I? Do we know that

> conditional simulation (permutation bootstrap) is valid in some cases, if

> so which? Is there a development of parametric bootstrap?

>

>

> Roger

>

>

>> best,

>>

>> Rafael HM Pereira

>> http://urbandemographics.blogspot.com

>>

>>

>> On Mon, Jul 31, 2017 at 10:52 PM, Roger Bivand <[hidden email]>

>> wrote:

>>

>> Rafael,

>>>

>>> I'm sorry, but there is no way you can logically "analyze who benefits

>>> the

>>> recent changes in the transport system in terms of access to jobs" from

>>> the

>>> data you have.

>>>

>>> Even if you had aggregate household income data for 2014 and 2017 (not

>>> for

>>> 2010 only), you would not know whether wealthier families had not

>>> dispaced

>>> poorer families as accessibility improved. You need individual data,

>>> either

>>> survey or register, preferably panel, to show that changes in

>>> accessibility

>>> change the economic welfare of households controlling for movement of

>>> households. The timestamps on the data make any attempt to do this very

>>> risky; the real findings from a hypothetical surevey-based panel might be

>>> completely different, especially if poorer households were displaced

>>> (also

>>> indirectly, through rising house prices driven by improved

>>> accessibility).

>>> Gauging the welfare effects of transport investments is very hard to

>>> instrument.

>>>

>>> The closest I could get was to map deciles of the change in access (more

>>> negatives than positives) and compare the aspatial income distributions:

>>>

>>> library(spdep)

>>> library(rgdal)

>>> map <- readOGR(dsn=".", layer="test_map")

>>> library(classInt)

>>> cI <- classIntervals(map$diffaccess, n=10, style="quantile")

>>> library(RColorBrewer)

>>> ybrpal <- brewer.pal(6, "YlOrBr")

>>> fC <- findColours(cI, ybrpal)

>>> qtm(map, fill="diffaccess", fill.breaks=cI$brks, format="Europe2")

>>> map$faccess <- factor(findInterval(map$diffaccess, cI$brks,

>>> all.inside=TRUE), labels=names(attr(fC, "table")))

>>> qtm(map, fill="diffaccess", fill.breaks=cI$brks, format="Europe2")

>>> acc_income <- split(map$income, map$faccess)

>>> do.call("rbind", lapply(acc_income, summary))

>>> dens <- lapply(acc_income, density)

>>> plot(1, ylab="", xlab="", type="n", xlim=c(-2000, 11000), ylim=c(0,

>>> 0.002))

>>> for (i in seq(along=dens)) lines(dens[[i]], col=i)

>>> legend("topright", legend=names(dens), col=1:length(dens), lty=1,

>>> bty="n")

>>>

>>> These density curves really do not suggest any clear relationship, other

>>> than that some areas with increased accessibility had higher incomes in

>>> 2010.

>>>

>>> You can examine the reverse relationship - were aggregate areas that were

>>> more wealthy in 2010 able to attract more changes to accessibility? The

>>> answer seems to be yes, they were able to do this:

>>>

>>> nb <- poly2nb(map)

>>> lw <- nb2listw(nb, style = "W", zero.policy = T)

>>> lm.morantest(lm(diffaccess ~ I(income/1000), map), lw)

>>> # SLX model

>>> summary(lmSLX(diffaccess ~ I(income/1000), map, lw))

>>> lm.morantest(lmSLX(diffaccess ~ I(income/1000), map, lw), lw)

>>> # Spatial Durbin error model - SDEM

>>> obj <- errorsarlm(diffaccess ~ I(income/1000), map, lw, etype="emixed")

>>> summary(impacts(obj))

>>> summary(impacts(lmSLX(diffaccess ~ I(income/1000), map, lw)))

>>> LR.sarlm(lmSLX(diffaccess ~ I(income/1000), map, lw), obj)

>>>

>>> It would be possible to run lm.morantest.sad() on the output of the SDEM

>>> model taking global spatial autocorrelation into account. If you need

>>> that,

>>> follow up in this thread.

>>>

>>> Bivariate Moran's I should not be used in this case, but could be used in

>>> other cases - use in change over time is troubling because randomisation

>>> will not be a good guide as t=1 and t=2 are subject to temporal as well

>>> as

>>> spatial autocorrelation, so you cannot use permutation bootstrap to find

>>> a

>>> usable measure of significance.

>>>

>>> Hope this clarifies, and thanks for the code.

>>>

>>> Roger

>>>

>>> On Sun, 30 Jul 2017, Rafael Pereira wrote:

>>>

>>> Roger,

>>>

>>>>

>>>> Population and income data are single point in time and come from the

>>>> 2010

>>>> Census.

>>>>

>>>> Accessibility variables in 2014 and 2017 show the proportion of jobs

>>>> accessible by public transport under 60 minutes. The variable diffaccess

>>>> shows the difference between these two. It's in percentage points

>>>> (access2017 - access2014)

>>>>

>>>> best,

>>>>

>>>> Rafael H M Pereira

>>>> urbandemographics.blogspot.com

>>>>

>>>> On Sun, Jul 30, 2017 at 7:41 AM, Roger Bivand <[hidden email]>

>>>> wrote:

>>>>

>>>> Thanks, I'll get back when able, offline now. What are the units of

>>>>

>>>>> observation, and are aggregate household incomes observed only once?

>>>>>

>>>>> Roger

>>>>>

>>>>> Roger Bivand

>>>>> Norwegian School of Economics

>>>>> Bergen, Norway

>>>>>

>>>>>

>>>>>

>>>>> Fra: Rafael Pereira

>>>>> Sendt: søndag 30. juli, 00.39

>>>>> Emne: Re: [R-sig-Geo] bivariate spatial correlation in R

>>>>> Kopi: Rogério Barbosa, [hidden email]

>>>>>

>>>>>

>>>>> Hi all, here is a reproducible example to calculate in R bivariate

>>>>> Moran's

>>>>> I and LISA clusters. This example is based on a this answer provided in

>>>>> SO*

>>>>> and it uses a toy model of my data. The R script and the shape file

>>>>> with

>>>>> the data are available on this link. https://gist.github.com/

>>>>> rafapereirabr/5348193abf779625f5e8c5090776a228 What this example does

>>>>> is

>>>>> to estimate the spatial association between household income per capita

>>>>> and

>>>>> the gains in accessibility to jobs. The aim is to analyze who benefits

>>>>> the

>>>>> recent changes in the transport system in terms of access to jobs. So

>>>>> the

>>>>> idea is not to find causal relationships, but spatial association

>>>>> between

>>>>> areas of high/low income who had high/low gains in accessibility. The

>>>>> variables in the data show info on the proportion of jobs accessible in

>>>>> both years 2014 and 2017 (access2014, access2017) and the difference

>>>>> between the two years in percentage points (diffaccess). Roger, I know

>>>>> you

>>>>> have shown to be a bit sceptical about this application of bivariate

>>>>> Moran's I. Do you still think a spatial regression would be more

>>>>> appropriate? Also, I would be glad to hear if others have comments on

>>>>> the

>>>>> code. This function is not implemented in any package so it would be

>>>>> great

>>>>> to have some feedback. Rafael H M Pereira

>>>>> urbandemographics.blogspot.com

>>>>> * https://stackoverflow.com/questions/45177590/map-of-

>>>>> bivariate-spatial-correlation-in-r-bivariate-lisa On Wed, Jul 26, 2017

>>>>> at

>>>>> 11:07 AM, Roger Bivand wrote: > On Wed, 26 Jul 2017, Rafael Pereira

>>>>> wrote:

>>>>>

>>>>> Roger, >> >> This example was provided only for the sake or making the

>>>>>>

>>>>>>>

>>>>>>> code easily >> reproducible for others and I'm more interested in how

>>>>>>

>>>>> the

>>>>> bi-variate >> Moran >> could be implemented in R, but your comments are

>>>>> very much welcomed and >> I've made changes to the question. >> >> My

>>>>> actual case study looks at bi-variate spatial correlation between (a)

>>>>> >>

>>>>> average household income per capita and (b) proportion of jobs in the

>>>>> city

>>>>>

>>>>> that are accessible under 60 minutes by transit. I don't think I could

>>>>>>

>>>>>>>

>>>>>>> use >> rates in this case but I will normalize the variables using >>

>>>>>>

>>>>> scale(data$variable). >> > > Please provide a reproducible example,

>>>>> either

>>>>> with a link to a data > subset, or using a builtin data set. My guess

>>>>> is

>>>>> that you do not need > bi-variate spatial correlation at all, but

>>>>> rather

>>>>> a

>>>>> spatial regression. > > The "causal" variable would then the the

>>>>> proportion

>>>>> of jobs accessible > within 60 minutes by transit, though this is

>>>>> extremely

>>>>> blunt, and lots of > other covariates (demography, etc.) impact average

>>>>> household income per > capita (per block/tract?). Since there are many

>>>>> missing variables in your > specification, any spatial correlation

>>>>> would

>>>>> be

>>>>> most closely associated > with them (demography, housing costs,

>>>>> education,

>>>>> etc.), and the choice of > units of measurement would dominate the

>>>>> outcome.

>>>>>

>>>>> This is also why bi-variate spatial correlation is seldom a good idea,

>>>>>>

>>>>>>>

>>>>>>> I > believe. It can be done, but most likely shouldn't, unless it can

>>>>>>

>>>>> be >

>>>>> motivated properly. > > By the way, the weighted and FDR-corrected SAD

>>>>> local Moran's I p-values of > the black/white ratio for Oregon (your

>>>>> toy

>>>>> example) did deliver the goods - > if you zoom in in mapview::mapview,

>>>>> you

>>>>> can see that it detects a rate > hotspot between the rivers. > > Roger

>>>>> >

>>>>>

>>>>>>

>>>>>>

>>>>> best, >> >> Rafael H M Pereira >> >> On Mon, Jul 24, 2017 at 7:56 PM,

>>>>>>

>>>>>>>

>>>>>>>> Roger Bivand >> wrote: >> >> On Mon, 24 Jul 2017, Rafael Pereira

>>>>>>>

>>>>>> wrote: >>>

>>>>>

>>>>> Hi all, >>> >>>> >>>> I would like to ask whether some you conducted

>>>>>>

>>>>>>>

>>>>>>>> bi-variate spatial >>>> correlation in R. >>>> >>>> I know the

>>>>>>>

>>>>>> bi-variate

>>>>> Moran's I is not implemented in the spdep library. >>>> I left a

>>>>> question

>>>>> on SO but also wanted to hear if anyone if the >>>> mainlist >>>> have

>>>>> come

>>>>> across this. >>>> https://stackoverflow.com/questions/45177590/map-of-

>>>>> bivariat >>>> e-spatial-correlation-in-r-bivariate-lisa >>>> >>>> I

>>>>> also

>>>>> know Roger Bivand has implemented the L index proposed by Lee >>>>

>>>>> (2001)

>>>>>

>>>>> in spdep, but I'm not I'm not sure whether the L local correlation

>>>>>>

>>>>>>> coefficients can be interpreted the same way as the local Moran's I

>>>>>>>>> coefficients. I couldn't find any reference commenting on this

>>>>>>>>> issue.

>>>>>>>>>

>>>>>>>>> I >>>> would very much appreciate your thoughts this. >>>> >>>> >>>

>>>>>>>>

>>>>>>> In the

>>>>> SO question, and in the follow-up, your presumably throw-away >>>

>>>>> example

>>>>> makes fundamental mistakes. The code in spdep by Virgilio >>>

>>>>> Gómez-Rubio

>>>>> is for uni- and bivariate L, and produces point values of >>> local >>>

>>>>> L.

>>>>> This isn't the main problem, which is rather that you are not taking

>>>>> >>>

>>>>> account of the underlying population counts, nor shrinking any

>>>>> estimates

>>>>>

>>>>> of >>> significance to accommodate population sizes. Population sizes

>>>>>>

>>>>>>>

>>>>>>>> vary from >>> 0 >>> to 11858, with the lower quartile at 3164 and

>>>>>>> upper

>>>>>>>

>>>>>> 5698: >>> plot(ecdf(oregon.tract$pop2000)). Should you be comparing

>>>>> rates

>>>>> in >>> stead? >>> These are also compositional variables (sum to

>>>>> pop2000,

>>>>> or 1 if rates) >>> with >>> the other missing components. You would

>>>>> probably be better served by >>> tools >>> examining spatial

>>>>> segregation,

>>>>> such as for example the seg package. >>> >>> The 0 count populations

>>>>> cause

>>>>> problems for an unofficial alternative, the >>> black/white ratio: >>>

>>>>>

>>>>>>

>>>>>>>> oregon.tract1 0,] >>> oregon.tract1$rat >> nb >> lw >> >>> which

>>>>> should

>>>>> still be adjusted by weighting: >>> >>> lm0 >> >>> I'm not advising

>>>>> this,

>>>>> but running localmoran.sad on this model output >>> yields SAD

>>>>> p-values <

>>>>> 0.05 after FDR correction only in contiguous tracts >>> on the

>>>>> Washington

>>>>> state line in Portland between the Columbia and >>> Williamette rivers.

>>>>> So

>>>>> do look at the variables you are using before >>> rushing into things.

>>>>>

>>>>>>

>>>>>>>>

>>>>> Hope this clarifies, >>> >>> Roger >>> >>> >>> best, >>>> >>>> Rafael

>>>>>>

>>>>>>>

>>>>>>>> HM Pereira >>>> http://urbandemographics.blogspot.com >>>> >>>>

>>>>>>>

>>>>>> [[alternative HTML version deleted]] >>>> >>>>

>>>>> _______________________________________________ >>>> R-sig-Geo mailing

>>>>> list >>>> [hidden email] >>>> 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 <+47%2055%2095%2093%2055>; e-mail:

>>>>> [hidden email] >>> Editor-in-Chief of The R Journal,

>>>>> https://journal.r-project.org/ >>> index.html >>>

>>>>> http://orcid.org/0000-0003-2392-6140 >>> https://scholar.google.no/

>>>>> citations?user=AWeghB0AAAAJ&hl=en >>> >> >> [[alternative HTML version

>>>>> deleted]] >> >> _______________________________________________ >>

>>>>> R-sig-Geo mailing list >> [hidden email] >>

>>>>> 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

>>>>> <+47%2055%2095%2093%2055>;

>>>>> e-mail: [hidden email] > Editor-in-Chief of The R Journal,

>>>>> https://journal.r-project.org/index.html > http://orcid.org/0000-0003-

>>>>> 2392-6140 > https://scholar.google.no/cita

>>>>> tions?user=AWeghB0AAAAJ&hl=en

>>>>>

>>>>>>

>>>>>> [[alternative HTML version deleted]] ______________________________

>>>>> _________________

>>>>> R-sig-Geo mailing list [hidden email]

>>>>> 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; e-mail: [hidden email]

>>> Editor-in-Chief of The R Journal, https://journal.r-project.org/

>>> index.html

>>> http://orcid.org/0000-0003-2392-6140

>>> https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

>>>

>>>

>>

> --

> Roger Bivand

> Department of Economics, Norwegian School of Economics,

> Helleveien 30, N-5045 Bergen, Norway.

> voice: +47 55 95 93 55; e-mail: [hidden email]

> Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html

> http://orcid.org/0000-0003-2392-6140

> https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

>

[[alternative HTML version deleted]]

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: How to draw the same legend (one legend) for the multiple spatial figures?

multiple spatial figures?

I would put the legend next to the upper-most plot and

wrap the whole thing in a single

/figure

so that it would be more obvious that the same legend applied to all.

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: writeOGR layer options for the KML driver?

> Hi,

>

> Where can I find documentation for available layer_options (if any) for

> the writeOGR KML driver?

http://www.gdal.org/drv_kml.html

for the KML driver,

http://www.gdal.org/drv_libkml.html

for the LIBKML driver. I an uncertain about their ability to do what you

need.

Roger

>

> Specifically, I want to be able to specify the line color when writing a

> SpatialLines or SpatialLinesDataFrame object to a KML file for viewing

> in Google Earth. (also specifyng a fill color and fill transparency

> would be nice, but I can certainly do without it)

>

> I see that the plotKML package has a kml() function. I'll send a

> separate email asking for assistance with it.

>

> Thanks

> -Don

>

>

--

Roger Bivand

Department of Economics, Norwegian School of Economics,

Helleveien 30, N-5045 Bergen, Norway.

voice: +47 55 95 93 55; e-mail: [hidden email]

Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html

http://orcid.org/0000-0003-2392-6140

https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Roger Bivand

Department of Economics

Norwegian School of Economics

Helleveien 30

N-5045 Bergen, Norway

### Help with plotKML::kml() ?

plotKML::kml()

function? I'm looking at various documentation and having trouble finding what I need.

Suppose "foo" is a SpatialPolygonsDataFrame object with a lat/long coordinate reference system. Then in R I can do, for example,

plot(foo, col='yellow', border='red')

I want to write a KML file that will display the polygon in Google Earth with the same colors.

plotKML::kml(foo, colour='yellow')

gets the fill, but I can't find how to set the border color.

I'd also like to be able to specify no fill, border only.

Thanks

-Don

--

Don MacQueen

Lawrence Livermore National Laboratory

7000 East Ave., L-627

Livermore, CA 94550

925-423-1062

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### writeOGR layer options for the KML driver?

Where can I find documentation for available layer_options (if any) for the writeOGR KML driver?

Specifically, I want to be able to specify the line color when writing a SpatialLines or SpatialLinesDataFrame object to a KML file for viewing in Google Earth.

(also specifyng a fill color and fill transparency would be nice, but I can certainly do without it)

I see that the plotKML package has a kml() function. I'll send a separate email asking for assistance with it.

Thanks

-Don

--

Don MacQueen

Lawrence Livermore National Laboratory

7000 East Ave., L-627

Livermore, CA 94550

925-423-1062

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### R-help: How to draw the same legend (one legend) for the multiple spatial figures?

# draw figure

## Open a new default device.

#get( getOption( "device" ) )()

## Set up plotting in two rows and three columns, plotting along rows first.

#par( mfrow = c( 3, 3 ) )

my.palette <- brewer.pal(n = 7, name = "OrRd") # set the color styles

## The first plot is located in row 1, column 1:

a<-spplot(nc1,sp.layout = list("sp.lines", as(lev, "SpatialLines")),main="All",col.regions=my.palette, cuts = 6, col = "transparent",scales=list(draw = TRUE),par.settings=list(fontsize=list(text=17)))

#plot( rnorm( n = 10 ), col = "red", main = "plot 1", cex.lab = 1.1 )

## The second plot is located in row 1, column 2:

b<-spplot(nc2,sp.layout = list("sp.lines", as(lev, "SpatialLines")),main="CC",col.regions=my.palette, cuts = 6, col = "transparent",scales=list(draw = TRUE),par.settings=list(fontsize=list(text=17)))

## The third plot is located in row 1, column 3:

c<-spplot(nc3,sp.layout = list("sp.lines", as(lev, "SpatialLines")),main="GF",col.regions=my.palette, cuts = 6, col = "transparent",scales=list(draw = TRUE),par.settings=list(fontsize=list(text=17)))

## The fourth plot is located in row 2, column 1:

d<-spplot(nc4,sp.layout = list("sp.lines", as(lev, "SpatialLines")),main="HA",col.regions=my.palette, cuts = 6, col = "transparent",scales=list(draw = TRUE),par.settings=list(fontsize=list(text=17)))

## The fifth plot is located in row 2, column 2:

e<-spplot(nc5,sp.layout = list("sp.lines", as(lev, "SpatialLines")),main="MI",col.regions=my.palette, cuts = 6, col = "transparent",scales=list(draw = TRUE),par.settings=list(fontsize=list(text=17)))

## The sixth plot is located in row 2, column 3:

f<-spplot(nc6,sp.layout = list("sp.lines", as(lev, "SpatialLines")),main="MP",col.regions=my.palette, cuts = 6, col = "transparent",scales=list(draw = TRUE),par.settings=list(fontsize=list(text=17)))

## The seventh plot is located in row 3, column 1:

g<-spplot(nc7,sp.layout = list("sp.lines", as(china, "SpatialLines")),main="MR",col.regions=my.palette, cuts = 6, col = "transparent",scales=list(draw = TRUE),par.settings=list(fontsize=list(text=17)))

Of them, nc1-nc7 are the raster data.

Thanks a lot!

Wei

[hidden email]

[[alternative HTML version deleted]]

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

### Re: Problem with localG in spdep

Roger Bivand

Norwegian School of Economics

Bergen, Norway

Fra: Victor Kuperman via R-sig-Geo

Sendt: torsdag 3. august, 23.26

Emne: [R-sig-Geo] Problem with localG in spdep

Til: [hidden email]

Dear all, I am new to using R for analyzing geographic data, so any input would be appreciated. My goal is to obtain the G* (local G) Getis-Ord statistic to a metric that is calculated for every US state and Canadian province (except islands, i.e. Hawaii and Prince Edward Island). The metric is standardized and its value range is [-2, 2]. Using spdep, I created a neighbors list from the US + Canada polygon data, and then generated a binary weight matrix. state_nb

[[alternative HTML version deleted]]

_______________________________________________

R-sig-Geo mailing list

[hidden email]

https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Roger Bivand

Department of Economics

Norwegian School of Economics

Helleveien 30

N-5045 Bergen, Norway