[localmoran.exact]

classic Classic list List threaded Threaded
3 messages Options
LM
Reply | Threaded
Open this post in threaded view
|

[localmoran.exact]

LM

Dear all,

I have been trying to obtain local Moran’s I with the exact approach (localmoran.exact)
 but I have some issues and would need your help.

I am trying to bring to light how spatial patterns of wealth have evolved from 2000 to 2016 in the European Union.

Looking at the Regional Gross Domestic Product (GDP) of European regions, my lm model looks like this :

lm(log(EU_NUTS@data[,year[i]])~1+factor(EU_NUTS@data$NUTS_L1))

where log(EU_NUTS@data[,year[i]]) are log of regional GDPs

and factor(EU_NUTS@data$NUTS_L1) is the vector of country dummies.

Sadly, I obtain an error message that says :

Error in integrate(integrand, lower = 0, upper = upper) :
  a limit is missing
Moreover : Warning messages:
1: In sqrt(2 * t2 - t1^2) : production of NaN
2: In sqrt(2 * t2 - t1^2) : production of NaN

According to the GitHub source code ,

        Vi <- listw2star <https://rdrr.io/cran/spdep/man/localmoran.sad.html>(B, select[i], style=style, n, D <https://rdrr.io/r/stats/deriv.html>, a,
            zero.policy=zero.policy)
        Viu <- lag.listw <https://rdrr.io/cran/spdep/man/lag.listw.html>(Vi, u, zero.policy=TRUE <https://rdrr.io/r/base/logical.html>)
        Ii <- c <https://rdrr.io/r/base/c.html>(crossprod <https://rdrr.io/r/base/crossprod.html>(u, Viu) / utu)
        ViX <- lag.listw <https://rdrr.io/cran/spdep/man/lag.listw.html>(Vi, X, zero.policy=TRUE <https://rdrr.io/r/base/logical.html>)
        MViM <- t <https://rdrr.io/r/base/t.html>(X) %*% ViX %*% XtXinv
        t1 <- -sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(MViM))
        sumsq.Vi <- function <https://rdrr.io/r/base/function.html>(x) {
            if <https://rdrr.io/r/base/Control.html> (is.null <https://rdrr.io/r/base/NULL.html>(x)) NA <https://rdrr.io/r/base/NA.html>
            else <https://rdrr.io/r/base/Control.html> sum <https://rdrr.io/r/base/sum.html>(x^2)
        }
        trVi2 <- sum <https://rdrr.io/r/base/sum.html>(sapply <https://rdrr.io/r/base/lapply.html>(Vi$weights <https://rdrr.io/r/stats/weights.html>, sumsq.Vi), na.rm=TRUE <https://rdrr.io/r/base/logical.html>)
        t2a <- sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(t <https://rdrr.io/r/base/t.html>(ViX) %*% ViX %*% XtXinv))
        t2b <- sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(MViM %*% MViM))
        t2 <- trVi2 - 2*t2a + t2b
        e1 <- 0.5 * (t1 + sqrt <https://rdrr.io/r/base/MathFun.html>(2*t2 - t1^2))
        en <- 0.5 * (t1 - sqrt <https://rdrr.io/r/base/MathFun.html>(2*t2 - t1^2))

My guess is that my setting causes (2*t2<t1^2) : Is it ? Why would it ?

My database is available here : http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt>
and my code : http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt>

I thank you very much for your attention and your help.

Best Regards
Lisa

PhD Student
GREDEG, Université Côté d’Azur






        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: [localmoran.exact]

Roger Bivand
Administrator
On Mon, 13 Aug 2018, Lisa Menez wrote:

>
> Dear all,
>
> I have been trying to obtain local Moran’s I with the exact approach
> (localmoran.exact) but I have some issues and would need your help.

A complete reproducible example is required (script and data, data may be
downloaded, do not attach). Does the same thing happen with
localmoran.sad()?

>
> I am trying to bring to light how spatial patterns of wealth have
> evolved from 2000 to 2016 in the European Union.
>
> Looking at the Regional Gross Domestic Product (GDP) of European regions, my lm model looks like this :
>
> lm(log(EU_NUTS@data[,year[i]])~1+factor(EU_NUTS@data$NUTS_L1))
>
> where log(EU_NUTS@data[,year[i]]) are log of regional GDPs
>
Never, ever use @ to access data. EU_NUTS@data$NUTS_L1 should be
EU_NUTS$NUTS_L1 and should be (made a) factor first. It isn't clear what
log(EU_NUTS@data[,year[i]]) will end up being - is "year" a character
vector? If so, EU_NUTS[[year[i]]] is the correct syntax.

> and factor(EU_NUTS@data$NUTS_L1) is the vector of country dummies.
>
> Sadly, I obtain an error message that says :
>
> Error in integrate(integrand, lower = 0, upper = upper) :
>  a limit is missing
> Moreover : Warning messages:
> 1: In sqrt(2 * t2 - t1^2) : production of NaN
> 2: In sqrt(2 * t2 - t1^2) : production of NaN
>
> According to the GitHub source code ,
No, as you can see from https://cran.r-project.org/package=spdep, the
development site of spdep is https://github.com/r-spatial/spdep/. Never
use non-authorised copies. To see the source, simply type the function
name. Use debug(localmoran.exact) to step through the function while it
runs, to see the values of the objects you are in doubt about. Do not do
this in RStudio, its debugging interface tries to be far too clever.

The following is illegible in plain text; post in plain text only.

Roger

>
>        Vi <- listw2star <https://rdrr.io/cran/spdep/man/localmoran.sad.html>(B, select[i], style=style, n, D <https://rdrr.io/r/stats/deriv.html>, a,
>    zero.policy=zero.policy)
>        Viu <- lag.listw <https://rdrr.io/cran/spdep/man/lag.listw.html>(Vi, u, zero.policy=TRUE <https://rdrr.io/r/base/logical.html>)
> Ii <- c <https://rdrr.io/r/base/c.html>(crossprod <https://rdrr.io/r/base/crossprod.html>(u, Viu) / utu)
>        ViX <- lag.listw <https://rdrr.io/cran/spdep/man/lag.listw.html>(Vi, X, zero.policy=TRUE <https://rdrr.io/r/base/logical.html>)
>        MViM <- t <https://rdrr.io/r/base/t.html>(X) %*% ViX %*% XtXinv
>        t1 <- -sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(MViM))
>        sumsq.Vi <- function <https://rdrr.io/r/base/function.html>(x) {
>            if <https://rdrr.io/r/base/Control.html> (is.null <https://rdrr.io/r/base/NULL.html>(x)) NA <https://rdrr.io/r/base/NA.html>
>    else <https://rdrr.io/r/base/Control.html> sum <https://rdrr.io/r/base/sum.html>(x^2)
> }
> trVi2 <- sum <https://rdrr.io/r/base/sum.html>(sapply <https://rdrr.io/r/base/lapply.html>(Vi$weights <https://rdrr.io/r/stats/weights.html>, sumsq.Vi), na.rm=TRUE <https://rdrr.io/r/base/logical.html>)
> t2a <- sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(t <https://rdrr.io/r/base/t.html>(ViX) %*% ViX %*% XtXinv))
> t2b <- sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(MViM %*% MViM))
> t2 <- trVi2 - 2*t2a + t2b
> e1 <- 0.5 * (t1 + sqrt <https://rdrr.io/r/base/MathFun.html>(2*t2 - t1^2))
> en <- 0.5 * (t1 - sqrt <https://rdrr.io/r/base/MathFun.html>(2*t2 - t1^2))
>
> My guess is that my setting causes (2*t2<t1^2) : Is it ? Why would it ?
>
> My database is available here : http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt>
> and my code : http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt>
>
> I thank you very much for your attention and your help.
>
> Best Regards
> Lisa
>
> PhD Student
> GREDEG, Université Côté d’Azur
>
>
>
>
>
>
> [[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]
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
LM
Reply | Threaded
Open this post in threaded view
|

Re: [localmoran.exact]

LM

> Le 14 août 2018 à 20:09, Roger Bivand <[hidden email]> a écrit :
>
> If you post, you then reply in-thread on the list. Otherwise, nobody else knows what the resolution was. :
Sorry. My bad.
>
> Your explanation below contains lots of inaccuracies - all Spatial*DataFrames have standard data.frame access methods as I showed.
No doubt about it ! I will try do better now that I know about it.  

Thank you Sir.
Lisa

>
> Roger Bivand
>
> On Tue, 14 Aug 2018, Lisa Menez wrote:
>
>> Sir Bivand,
>>
>> The complete example was here :
>>
>> My database is available here : http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt> and my code : http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt>
>>
>> The localmoran.sad() actually works. Is it because « If we assume global spatial independence, we can give for any projection matrix M the spectrum of eigenvalues as well as the saddle point approximation of local Moran’s Ii in analytical terms. This increases the computational efficiency as it avoids reverting to numerical calculations of the spectrum of eigenvalues and finding iteratively the root of the saddle point equation » ?
>>
>> To the best of my knowledge, I am constrained by the settings of my data (extracted from Eurostat) on calling EU_NUTS@data rather than EU_NUTS. Indeed EU_NUTS are spatial polygons and EU_NUTS@data their relative attributes. I failed to change this name.
>>
>> for each i, EU_NUTS@data[,year[i]] is a numerical vector. year is a character vector containing the sample years, allowing to go through the corresponding EU_NUTS@data' columns.
>>
>> On your advice, I tried to go through debug(localmoran.exact). As I could not find any good reason why (2*t2<t1^2), I turned to your last remark : Changing NUTS_L1 into a factor before using it into the loop. It actually was the issue. I miscoded the factor, keeping character chains and not transforming them into numbers.
>>
>>
>> Thank you very much for your help and recommandations.
>> Lisa
>>
>>
>>> Le 13 août 2018 à 16:12, Roger Bivand <[hidden email]> a écrit :
>>>
>>> On Mon, 13 Aug 2018, Lisa Menez wrote:
>>>
>>>>
>>>> Dear all,
>>>>
>>>> I have been trying to obtain local Moran’s I with the exact approach (localmoran.exact) but I have some issues and would need your help.
>>>
>>> A complete reproducible example is required (script and data, data may be downloaded, do not attach). Does the same thing happen with localmoran.sad()?
>>>
>>>>
>>>> I am trying to bring to light how spatial patterns of wealth have evolved from 2000 to 2016 in the European Union.
>>>>
>>>> Looking at the Regional Gross Domestic Product (GDP) of European regions, my lm model looks like this :
>>>>
>>>> lm(log(EU_NUTS@data[,year[i]])~1+factor(EU_NUTS@data$NUTS_L1))
>>>>
>>>> where log(EU_NUTS@data[,year[i]]) are log of regional GDPs
>>>>
>>>
>>> Never, ever use @ to access data. EU_NUTS@data$NUTS_L1 should be EU_NUTS$NUTS_L1 and should be (made a) factor first. It isn't clear what log(EU_NUTS@data[,year[i]]) will end up being - is "year" a character vector? If so, EU_NUTS[[year[i]]] is the correct syntax.
>>>
>>>> and factor(EU_NUTS@data$NUTS_L1) is the vector of country dummies.
>>>>
>>>> Sadly, I obtain an error message that says :
>>>>
>>>> Error in integrate(integrand, lower = 0, upper = upper) :
>>>> a limit is missing
>>>> Moreover : Warning messages:
>>>> 1: In sqrt(2 * t2 - t1^2) : production of NaN
>>>> 2: In sqrt(2 * t2 - t1^2) : production of NaN
>>>>
>>>> According to the GitHub source code ,
>>>
>>> No, as you can see from https://cran.r-project.org/package=spdep, the development site of spdep is https://github.com/r-spatial/spdep/. Never use non-authorised copies. To see the source, simply type the function name. Use debug(localmoran.exact) to step through the function while it runs, to see the values of the objects you are in doubt about. Do not do this in RStudio, its debugging interface tries to be far too clever.
>>>
>>> The following is illegible in plain text; post in plain text only.
>>>
>>> Roger
>>>
>>>>
>>>>      Vi <- listw2star <https://rdrr.io/cran/spdep/man/localmoran.sad.html>(B, select[i], style=style, n, D <https://rdrr.io/r/stats/deriv.html>, a,
>>>>    zero.policy=zero.policy)
>>>>      Viu <- lag.listw <https://rdrr.io/cran/spdep/man/lag.listw.html>(Vi, u, zero.policy=TRUE <https://rdrr.io/r/base/logical.html>)
>>>> Ii <- c <https://rdrr.io/r/base/c.html>(crossprod <https://rdrr.io/r/base/crossprod.html>(u, Viu) / utu)
>>>>      ViX <- lag.listw <https://rdrr.io/cran/spdep/man/lag.listw.html>(Vi, X, zero.policy=TRUE <https://rdrr.io/r/base/logical.html>)
>>>>      MViM <- t <https://rdrr.io/r/base/t.html>(X) %*% ViX %*% XtXinv
>>>>      t1 <- -sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(MViM))
>>>>      sumsq.Vi <- function <https://rdrr.io/r/base/function.html>(x) {
>>>>          if <https://rdrr.io/r/base/Control.html> (is.null <https://rdrr.io/r/base/NULL.html>(x)) NA <https://rdrr.io/r/base/NA.html>
>>>>    else <https://rdrr.io/r/base/Control.html> sum <https://rdrr.io/r/base/sum.html>(x^2)
>>>> }
>>>> trVi2 <- sum <https://rdrr.io/r/base/sum.html>(sapply <https://rdrr.io/r/base/lapply.html>(Vi$weights <https://rdrr.io/r/stats/weights.html>, sumsq.Vi), na.rm=TRUE <https://rdrr.io/r/base/logical.html>)
>>>> t2a <- sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(t <https://rdrr.io/r/base/t.html>(ViX) %*% ViX %*% XtXinv))
>>>> t2b <- sum <https://rdrr.io/r/base/sum.html>(diag <https://rdrr.io/r/base/diag.html>(MViM %*% MViM))
>>>> t2 <- trVi2 - 2*t2a + t2b
>>>> e1 <- 0.5 * (t1 + sqrt <https://rdrr.io/r/base/MathFun.html>(2*t2 - t1^2))
>>>> en <- 0.5 * (t1 - sqrt <https://rdrr.io/r/base/MathFun.html>(2*t2 - t1^2))
>>>>
>>>> My guess is that my setting causes (2*t2<t1^2) : Is it ? Why would it ?
>>>>
>>>> My database is available here : http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/mydata2018-1.txt>
>>>> and my code : http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt <http://hcc.cnrs.fr/wp-content/uploads/2018/08/myCode_2018.txt>
>>>>
>>>> I thank you very much for your attention and your help.
>>>>
>>>> Best Regards
>>>> Lisa
>>>>
>>>> PhD Student
>>>> GREDEG, Université Côté d’Azur
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> [[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]
>>> 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]
> 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