Dear colleagues,
I have a 500 by 500 matrix with presence or absence value in each cell as an input matrix. I am trying to create an output matrix having in each cell the sum of the values from neigbouring cells in the input matrix according to a moving-window defined as: 0 1 1 1 0 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 0 1 1 1 0 where the output cell receiving the values is located in the center. I have written a function that loop over the cells of the output matrix, and sum up the values from these neigbouring cells. It works, but it is fairly slow, and hangs if I try it over a larger matrix. Given that this is a rather basic function available in most raster GIS packages (IDRISI and AV Spatial Analyst) and running fairly quickly, I wondered whether anyone had previously come accross a similar function in an existing R-package. Any helps or pointer would be appreciated, Many thanks in advance, Marius |
Hello,
I think that you could do it using spdep package. I suppose that, if you: - construct a neighborhoord object (nb class). You can define a distance based neighborhood (using dnearneigh) as it seems that your definition is similar to say i and j neighbors if 0<dij<sqrt(cellsize^2+(2*cellsize)^2) - then transform the nb to listw with nb2listw with style="B" - then use lag.listw to sompute the sum an example: > myvalues=rbinom(100,size=1,prob=0.5) > xy=expand.grid(1:10,1:10) > dnn1=dnearneigh(as.matrix(xy),0,sqrt(1+4)) > plot(dnn1,as.matrix(xy)) > nb2listw(dnn1,style="B") > lw1=nb2listw(dnn1,style="B") > sum.neigh=lag.listw(lw1,myvalues) > sum.neigh [1] 3 6 4 5 7 5 4 4 3 3 6 8 6 7 9 8 6 8 6 4 5 9 10 11 10 7 8 9 7 5 [31] 7 10 14 11 10 8 8 9 8 5 8 12 11 10 11 12 10 12 11 7 9 11 12 9 13 11 10 12 12 9 [61] 6 11 12 13 10 12 10 13 9 7 6 9 12 10 12 10 13 12 9 5 6 8 8 10 10 11 8 9 7 5 [91] 4 6 7 6 9 7 6 4 5 3 Roger could confirm that it is correct ! At 10:05 16/12/2004, Marius Gilbert wrote: >Dear colleagues, > >I have a 500 by 500 matrix with presence or absence value in each cell as >an input matrix. I am trying to create an output matrix having in each >cell the sum of the values from neigbouring cells in the input matrix >according to a moving-window defined as: > >0 1 1 1 0 >1 1 1 1 1 >1 1 0 1 1 >1 1 1 1 1 >0 1 1 1 0 > >where the output cell receiving the values is located in the center. I >have written a function that loop over the cells of the output matrix, and >sum up the values from these neigbouring cells. It works, but it is fairly >slow, and hangs if I try it over a larger matrix. Given that this is a >rather basic function available in most raster GIS packages (IDRISI and AV >Spatial Analyst) and running fairly quickly, I wondered whether anyone had >previously come accross a similar function in an existing R-package. > >Any helps or pointer would be appreciated, > >Many thanks in advance, > >Marius > >_______________________________________________ >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/ |
Stephane, have you tried how well this scales to, say, 1000x1000
matrices? Indeed, most basic raster GIS functionality is not available in R packages, and would be welcome. R code is efficient for small problems, but for cases like this scales badly as Marius noted. We are actually working on an interface from R to the functionality of PCRaster (http://pcraster.geog.uu.nl -- yet another raster GIS), and expect to have this working within a few months. PCRaster is not open source, so it is unlikely to become a CRAN package straight away, but still it may fill a need. This all depends on availability of sp, the long awaited package with S classes for spatial point, line, grid, arc/node and polygon data. We're working very hard on this, and most likely will announce a beta in January. I expect a PCRasteR package to follow that soon. -- Edzer Stephane DRAY wrote: > Hello, > I think that you could do it using spdep package. I suppose that, if you: > - construct a neighborhoord object (nb class). You can define a distance > based neighborhood (using dnearneigh) as it seems that your definition > is similar to say i and j neighbors if > 0<dij<sqrt(cellsize^2+(2*cellsize)^2) > - then transform the nb to listw with nb2listw with style="B" > - then use lag.listw to sompute the sum > > an example: > > myvalues=rbinom(100,size=1,prob=0.5) > > xy=expand.grid(1:10,1:10) > > dnn1=dnearneigh(as.matrix(xy),0,sqrt(1+4)) > > plot(dnn1,as.matrix(xy)) > > nb2listw(dnn1,style="B") > > lw1=nb2listw(dnn1,style="B") > > sum.neigh=lag.listw(lw1,myvalues) > > sum.neigh > [1] 3 6 4 5 7 5 4 4 3 3 6 8 6 7 9 8 6 8 6 4 5 9 > 10 11 10 7 8 9 7 5 > [31] 7 10 14 11 10 8 8 9 8 5 8 12 11 10 11 12 10 12 11 7 9 11 > 12 9 13 11 10 12 12 9 > [61] 6 11 12 13 10 12 10 13 9 7 6 9 12 10 12 10 13 12 9 5 6 8 > 8 10 10 11 8 9 7 5 > [91] 4 6 7 6 9 7 6 4 5 3 > > Roger could confirm that it is correct ! > > > > At 10:05 16/12/2004, Marius Gilbert wrote: > >> Dear colleagues, >> >> I have a 500 by 500 matrix with presence or absence value in each cell as >> an input matrix. I am trying to create an output matrix having in each >> cell the sum of the values from neigbouring cells in the input matrix >> according to a moving-window defined as: >> >> 0 1 1 1 0 >> 1 1 1 1 1 >> 1 1 0 1 1 >> 1 1 1 1 1 >> 0 1 1 1 0 >> >> where the output cell receiving the values is located in the center. I >> have written a function that loop over the cells of the output matrix, >> and >> sum up the values from these neigbouring cells. It works, but it is >> fairly >> slow, and hangs if I try it over a larger matrix. Given that this is a >> rather basic function available in most raster GIS packages (IDRISI >> and AV >> Spatial Analyst) and running fairly quickly, I wondered whether anyone >> had >> previously come accross a similar function in an existing R-package. >> >> Any helps or pointer would be appreciated, >> >> Many thanks in advance, >> >> Marius >> >> _______________________________________________ >> 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/ > > _______________________________________________ > R-sig-Geo mailing list > R-sig-Geo at stat.math.ethz.ch > https://stat.ethz.ch/mailman/listinfo/r-sig-geo |
Does anyone know of a GIS based on a regular hexagon tessellation rather
than on squares or rectangles? (There was a project twenty some years ago called System X developed for the Navy by Dean Lucas that started into commercialization but never finished as far as I know.) r-sig-geo-bounces at stat.math.ethz.ch wrote on 2004-12-16 12:17:28: > Stephane, have you tried how well this scales to, say, 1000x1000 > matrices? Indeed, most basic raster GIS functionality is not available > in R packages, and would be welcome. R code is efficient for small > problems, but for cases like this scales badly as Marius noted. > > We are actually working on an interface from R to the functionality of > PCRaster (http://pcraster.geog.uu.nl -- yet another raster GIS), and > expect to have this working within a few months. PCRaster is not open > source, so it is unlikely to become a CRAN package straight away, but > still it may fill a need. > > This all depends on availability of sp, the long awaited package with S > classes for spatial point, line, grid, arc/node and polygon data. We're > working very hard on this, and most likely will announce a beta in > January. I expect a PCRasteR package to follow that soon. > -- > Edzer |
Free forum by Nabble | Edit this page |