[R-sig-Geo] problems with dnearneigh of spdep
Roger Bivand
Roger.Bivand at nhh.no
Sat Oct 30 12:53:16 CEST 2004
Stéphane:
A test new release source package is on:
http://spatial.nhh.no/R/spdep/spdep_0.2-23.tar.gz
and a test Windows binary package built for R 2.0.0:
http://spatial.nhh.no/R/spdep/spdep_0.2-23.zip
If you are on R < 2.0.0 and need to carry on using the Windows binary
package you have, please edit the dnearneigh() function as follows:
if (!lonlat) {
for (i in 1:dimension) md <- sum(md, (diff(range(x[,i]))^2))
if (d2 > sqrt(md)) d2 <- sqrt(md)
}
to:
if (!lonlat) {
for (i in 1:dimension) md <- sum(md, (diff(range(x[,i]))^2))
md <- md + sqrt(.Machine$double.eps)
if (d2 > sqrt(md)) d2 <- sqrt(md)
}
The additional line adds "fuzz" to the "sane" maximum distance to let
tests of dist <= d2 pass when they are equal - I think some truncation
beyond reasonable numerical precision is what is going on, and this sanity
test for unreasonable distances then goes wrong. I've tried up to 40x40
grids on i386 Linux and i386 Windows, now all numbers agree.
Please let me know if this fixed the problem for you,
Roger
On Sat, 30 Oct 2004, Roger Bivand wrote:
> On Fri, 29 Oct 2004, Roger Bivand wrote:
>
> > On Fri, 29 Oct 2004, Stephane DRAY wrote:
> >
> > > My version :
> > >
> > > > version
> > > _
> > > platform i386-pc-mingw32
> > > arch i386
> > > os mingw32
> > > system i386, mingw32
> > > status
> > > major 1
> > > minor 9.1
> > > year 2004
> > > month 06
> > > day 21
> > > language R
> > > > spdep()
> > > [1] "spdep, version 0.2-18, 2004-05-25"
> > > > sum(card(mydnn))
> > > [1] 9896
> > > > all(sapply(mydnn, length) == 99)
> > > [1] FALSE
> >
> > I can replicate on Windows, R 1.9.1, spdep 0.2-18. I'll report back on R
> > 2.0.0, spdep 0.2-22 in a little while.
>
> It also replicates with fresh binaries on:
>
> > spdep()
> [1] "spdep, version 0.2-22, 2004-08-07"
> > version
> _
> platform i386-pc-mingw32
> arch i386
> os mingw32
> system i386, mingw32
> status
> major 2
> minor 0.0
> year 2004
> month 10
> day 04
> language R
> > mydnn
> Neighbour list object:
> Number of regions: 100
> Number of nonzero links: 9896
> Percentage nonzero weights: 98.96
> Average number of links: 98.96
>
> I'll look and see what is going on, thanks for the bug report.
>
> Roger
>
> >
> > Roger
> >
> > >
> > > It is not the last one, I will try to update later. I think that the
> > > problem is in the calculation of neighbors. Perhaps, the fortran code has
> > > changed. I do not have the sources, so I could not verify.
> > >
> > > Thanks !
> > >
> > >
> > >
> > > At 13:01 29/10/2004, Roger Bivand wrote:
> > > >On Fri, 29 Oct 2004, Stephane DRAY wrote:
> > > >
> > > > > Hello,
> > > > > I have some strange results with the dnearneigh function.
> > > > > I would like to have all points of grid linked. I do something like that:
> > > > >
> > > > > > xy=expand.grid(1:10,1:10)
> > > > > > mydnn=dnearneigh(as.matrix(xy),d1=0,d2=50)
> > > > > > mydnn
> > > > > Neighbour list object:
> > > > > Number of regions: 100
> > > > > Number of nonzero links: 9896
> > > > > Percentage nonzero weights: 98.96
> > > > > Average number of links: 98.96
> > > > > I assume to have 100*100-100=9900 links.
> > > > >
> > > > > > m01=nb2mat(mydnn,style="B")
> > > > > > apply(m01,1,sum)
> > > > > 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
> > > > > 19 20 21 22 23 24 25 26 27
> > > > > 98 99 99 99 99 99 99 99 99 98 99 99 99 99 99 99 99 99
> > > > > 99 99 99 99 99 99 99 99 99
> > > > > 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45
> > > > > 46 47 48 49 50 51 52 53 54
> > > > > 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99
> > > > > 99 99 99 99 99 99 99 99 99
> > > > > 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72
> > > > > 73 74 75 76 77 78 79 80 81
> > > > > 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99
> > > > > 99 99 99 99 99 99 99 99 99
> > > > > 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98
> > > > 99 100
> > > > > 99 99 99 99 99 99 99 99 99 98 99 99 99 99 99 99 99
> > > > 99 98
> > > > >
> > > > > It seems that the first and the last points are not linked. I have
> > > > choose a
> > > > > value of d2 very high, and the distance between this two points is only
> > > > > sqrt(162)
> > > >
> > > >I can't replicate this:
> > > >
> > > > > version
> > > > _
> > > >platform i686-pc-linux-gnu
> > > >arch i686
> > > >os linux-gnu
> > > >system i686, linux-gnu
> > > >status
> > > >major 2
> > > >minor 0.0
> > > >year 2004
> > > >month 10
> > > >day 04
> > > >language R
> > > > > library(spdep)
> > > > > spdep()
> > > >[1] "spdep, version 0.2-22, 2004-08-07"
> > > > > xy=expand.grid(1:10,1:10)
> > > > > mydnn=dnearneigh(as.matrix(xy),d1=0,d2=50)
> > > > > mydnn
> > > >Neighbour list object:
> > > >Number of regions: 100
> > > >Number of nonzero links: 9900
> > > >Percentage nonzero weights: 99
> > > >Average number of links: 99
> > > >
> > > >The line "Number of nonzero links:" is calculated as:
> > > >
> > > > > sum(card(mydnn))
> > > >[1] 9900
> > > >
> > > >You can also look at:
> > > >
> > > > > all(sapply(mydnn, length) == 99)
> > > >[1] TRUE
> > > > > which(sapply(mydnn, length) != 99)
> > > >numeric(0)
> > > >
> > > >if you want to be sure it isn't card() getting it wrong. Puzzling. I don't
> > > >think any of the relevant code has changed since September 2003.
> > > >
> > > >Roger
> > > >
> > > > >
> > > > > Any ideas ?
> > > > > Thanks
> > > > > Stéphane DRAY
> > > > >
> > > > --------------------------------------------------------------------------------------------------
> > > >
> > > > >
> > > > > Département des Sciences Biologiques
> > > > > Université de Montréal, C.P. 6128, succursale centre-ville
> > > > > Montréal, Québec H3C 3J7, Canada
> > > > >
> > > > > Tel : (514) 343-6111 poste 1233 Fax : (514) 343-2293
> > > > > E-mail : stephane.dray at umontreal.ca
> > > > >
> > > > --------------------------------------------------------------------------------------------------
> > > >
> > > > >
> > > > >
> > > > Web http://www.steph280.freesurf.fr/
> > > > >
> > > > > _______________________________________________
> > > > > R-sig-Geo mailing list
> > > > > R-sig-Geo at stat.math.ethz.ch
> > > > > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> > > > >
> > > >
> > > >--
> > > >Roger Bivand
> > > >Economic Geography Section, Department of Economics, Norwegian School of
> > > >Economics and Business Administration, Breiviksveien 40, N-5045 Bergen,
> > > >Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93
> > > >e-mail: Roger.Bivand at nhh.no
> > > >
> > > >_______________________________________________
> > > >R-sig-Geo mailing list
> > > >R-sig-Geo at stat.math.ethz.ch
> > > >https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> > >
> > > Stéphane DRAY
> > > --------------------------------------------------------------------------------------------------
> > >
> > > Département des Sciences Biologiques
> > > Université de Montréal, C.P. 6128, succursale centre-ville
> > > Montréal, Québec H3C 3J7, Canada
> > >
> > > Tel : (514) 343-6111 poste 1233 Fax : (514) 343-2293
> > > E-mail : stephane.dray at umontreal.ca
> > > --------------------------------------------------------------------------------------------------
> > >
> > > Web http://www.steph280.freesurf.fr/
> > > --------------------------------------------------------------------------------------------------
> > >
> > >
> > >
> >
> >
>
>
--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Breiviksveien 40, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93
e-mail: Roger.Bivand at nhh.no
More information about the R-sig-Geo
mailing list