nearest distance in matrix

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

nearest distance in matrix

Marius Gilbert
Hello,

I'm trying to use R to simulate biological invasions, and got stuck with
the following:

Is there a function that uses a matrix of occupied/empty values (1 or 0)
as input, and producing an output matrix of the same size, with each cell
containing the distance to the nearest occupied cell (1) of the input
matrix ?

For example, the input matrix:

0    1    1    0
0    0    1    0
0    0    0    0
0    0    0    0

Would produce an output matrix like:

1    0    0    1
1.4  1    0    1
2.2  1.4  1    1.4
2.8  2.2  2.0  2.2

For those familiar with AV Spatial Analyst, this function would be the
equivalent of OutGrid = EucDistance (directionFN, allocationFN,
maxDistance).

Many thanks in advance for any help or pointers,

Marius



Reply | Threaded
Open this post in threaded view
|

nearest distance in matrix

Barry Rowlingson
Marius Gilbert wrote:
> Hello,
>
> I'm trying to use R to simulate biological invasions, and got stuck with
> the following:
>
> Is there a function that uses a matrix of occupied/empty values (1 or 0)
> as input, and producing an output matrix of the same size, with each cell
> containing the distance to the nearest occupied cell (1) of the input
> matrix ?

  Try this, which uses an internal, undocumented Spatstat function
'exactPdt':

library(spatstat)
nmatdist <- function(aMat){
   storage.mode(aMat) <- "logical"
   win <- owin(xrange=c(1,ncol(aMat)),yrange=c(1,nrow(aMat)),mask=aMat)
   m <- exactPdt(win)
   return(m$d)
}

> For example, the input matrix:
>
> 0    1    1    0
> 0    0    1    0
> 0    0    0    0
> 0    0    0    0
>

> Would produce an output matrix like:
>
> 1    0    0    1
> 1.4  1    0    1
> 2.2  1.4  1    1.4
> 2.8  2.2  2.0  2.2

  test:

 > m1
      [,1] [,2] [,3] [,4] [,5]
[1,]    0    0    0    0    0
[2,]    1    0    0    0    0
[3,]    1    1    0    0    0
[4,]    0    0    0    0    0
 > nmatdist(m1)
      [,1]     [,2]     [,3]     [,4]     [,5]
[1,]    1 1.414214 2.236068 2.828427 3.605551
[2,]    0 1.000000 1.414214 2.236068 3.162278
[3,]    0 0.000000 1.000000 2.000000 3.000000
[4,]    1 1.000000 1.414214 2.236068 3.162278

  looks about right!

Baz



Reply | Threaded
Open this post in threaded view
|

nearest distance in matrix

White.Denis@epamail.epa.gov
In reply to this post by Marius Gilbert




Here's another version that probably can be simplified:

nmatdist <- function (m)
{
    v <- as.vector (m)
    ind <- as.matrix (expand.grid (seq (nrow (m)), seq (ncol (m))))
    d <- as.matrix (dist (ind, upper=TRUE))
    ones <- as.numeric (dimnames(ind[v == 1,])[[1]])
    matrix (sapply (seq(nrow(ind)),
        function (x) min (d[x, ones])), nrow=nrow(m))
}

Denis White
   US EPA, 200 SW 35th St, Corvallis, Oregon, 97333 USA
   voice: 541.754.4476, email: white.denis at epa.gov
   web: www.epa.gov/wed/pages/staff/white/

r-sig-geo-bounces at stat.math.ethz.ch wrote on 2004-07-14 09:50:16:

> Hello,
>
> I'm trying to use R to simulate biological invasions, and got stuck
with
> the following:
>
> Is there a function that uses a matrix of occupied/empty values (1 or
0)
> as input, and producing an output matrix of the same size, with each
cell

> containing the distance to the nearest occupied cell (1) of the input
> matrix ?
>
> For example, the input matrix:
>
> 0    1    1    0
> 0    0    1    0
> 0    0    0    0
> 0    0    0    0
>
> Would produce an output matrix like:
>
> 1    0    0    1
> 1.4  1    0    1
> 2.2  1.4  1    1.4
> 2.8  2.2  2.0  2.2
>
> For those familiar with AV Spatial Analyst, this function would be the

> equivalent of OutGrid = EucDistance (directionFN, allocationFN,
> maxDistance).
>
> Many thanks in advance for any help or pointers,
>
> Marius
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://www.stat.math.ethz.ch/mailman/listinfo/r-sig-geo



Reply | Threaded
Open this post in threaded view
|

nearest distance in matrix

Barry Rowlingson
White.Denis at epamail.epa.gov wrote:

>
>
>
> Here's another version that probably can be simplified:
>
> nmatdist <- function (m)
> {
>     v <- as.vector (m)
>     ind <- as.matrix (expand.grid (seq (nrow (m)), seq (ncol (m))))
>     d <- as.matrix (dist (ind, upper=TRUE))
>     ones <- as.numeric (dimnames(ind[v == 1,])[[1]])
>     matrix (sapply (seq(nrow(ind)),
>         function (x) min (d[x, ones])), nrow=nrow(m))
> }

  Eeek! Have you tried that on a 100x100 matrix? Darn near killed my
machine!

  > mt=matrix(runif(10000)>.1,100,100)
  > nmatdist(mt)
   ... near swap death experience...

  The code used by Adrian Baddeley's spatstat routine uses a very neat
method for working out the distances, which involves sweeping along rows
and columns or something. He did explain it to me when I was in Perth
but I can't recall it now!

  Anyway, its super-quick and uses next-to-no memory. Here's how long my
function that calls the spatstat routine takes:

  > unix.time(nmatdist(mt))
  [1] 0.02 0.01 0.02 0.00 0.00

  it was so quick I thought I'd check I'd not done it on a small matrix
by mistake:

  > dim(mt)
  [1] 100 100

Nope!

Baz



Reply | Threaded
Open this post in threaded view
|

nearest distance in matrix

White.Denis@epamail.epa.gov




Yup.  Space complexity is O((nrow(mat)*ncol(mat))^2).


r-sig-geo-bounces at stat.math.ethz.ch wrote on 2004-07-14 10:49:14:

> White.Denis at epamail.epa.gov wrote:
> >
> >
> >
> > Here's another version that probably can be simplified:
> >
> > nmatdist <- function (m)
> > {
> >     v <- as.vector (m)
> >     ind <- as.matrix (expand.grid (seq (nrow (m)), seq (ncol (m))))
> >     d <- as.matrix (dist (ind, upper=TRUE))
> >     ones <- as.numeric (dimnames(ind[v == 1,])[[1]])
> >     matrix (sapply (seq(nrow(ind)),
> >         function (x) min (d[x, ones])), nrow=nrow(m))
> > }
>
>   Eeek! Have you tried that on a 100x100 matrix? Darn near killed my
> machine!
>
>   > mt=matrix(runif(10000)>.1,100,100)
>   > nmatdist(mt)
>    ... near swap death experience...
>
>   The code used by Adrian Baddeley's spatstat routine uses a very neat

> method for working out the distances, which involves sweeping along
rows
> and columns or something. He did explain it to me when I was in Perth
> but I can't recall it now!
>
>   Anyway, its super-quick and uses next-to-no memory. Here's how long
my
> function that calls the spatstat routine takes:
>
>   > unix.time(nmatdist(mt))
>   [1] 0.02 0.01 0.02 0.00 0.00
>
>   it was so quick I thought I'd check I'd not done it on a small
matrix

> by mistake:
>
>   > dim(mt)
>   [1] 100 100
>
> Nope!
>
> Baz
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://www.stat.math.ethz.ch/mailman/listinfo/r-sig-geo



Reply | Threaded
Open this post in threaded view
|

nearest distance in matrix

Marius Gilbert
In reply to this post by Barry Rowlingson
Dear Barry,

>   The code used by Adrian Baddeley's spatstat routine uses a very neat
> method for working out the distances, which involves sweeping along rows
> and columns or something. He did explain it to me when I was in Perth
> but I can't recall it now!

I knew there should be way to do it without estimating the full distance
redistribution matrix O((nrow(mat)*ncol(mat))^2) because the equivalent  
AV Spatial analyst can do rather big raster (3000 x 3000 pixels) in a few
seconds.

>   Anyway, its super-quick and uses next-to-no memory. Here's how long my
> function that calls the spatstat routine takes:

It is super quick indeed, I tested it on a 1600 x 250 matrix, ang got the
output 1600x250 distance matrix in less than 5 seconds on a PIV 1700 Mhz,
so fairly efficient indeed,

Many thanks, there's no I could have found this myself,

Cheers,

Marius