[postgis-users] spatial distribution maps (heat maps?)

[Phil James]

Sorry there was an error in that code it should be:

import rpy2.robjects as robjects
import sys
#Create R stuff and load libraries
R = robjects.rs
#R.library('sp')
R.library('gstat')
R.library('raster')
R.library('RPostgreSQL')
R("drv <- dbDriver('PostgreSQL')")
R("con <- dbConnect(drv, dbname='progsci', user='postgres')")
R("rs <- dbSendQuery(con, 'SELECT * FROM crimes WHERE crimetype='Burglary')")
R("data <- fetch(rs, n=-1)")
R("coordinates(data) <- ~x+y")
R("x.range <- as.integer(range(e at coords[,1]))")
R("y.range <- as.integer(range(e at coords[,2]))")
R("grd <- expand.grid(x=seq(from=x.range[1], to=x.range[2], bx=25), y=seq(from=y.range[1], to=y.range[2], by=25))")
R("coordinates(grd) <- ~x+y")
R("gridded(grd) <- TRUE")
R("idw.out <- idw(value~1, e, grd, idp=2.5)")
spplottext = "spplot(idw.out, \"var1.pred\", col.regions=topo.colors(20), pretty=TRUE)"
print(spplottext)
res = R(spplottext)
print(res)
sys.exit(0)



>________________________________
> From: Phil James <borntopedal at yahoo.co.uk>
>To: PostGIS Users Discussion <postgis-users at postgis.refractions.net> 
>Sent: Tuesday, 6 December 2011, 7:58
>Subject: Re: [postgis-users] spatial distribution maps (heat maps?)
> 
>
>I would use something that is designed to do this analysis - R would be an obvious choice but there are others - using Python you can connect to Postgis to grab the data and then rpy to run R commands from python - we have used this configuration successfully using Django and it works fast enough for the web generating an image - this could of course be optimised to cache outputs if performance is an issue.
>
>
>A synopsis of spatial functions in R is shown below:
>
>
>http://cran.r-project.org/web/views/Spatial.html
>
>
>
>Some sample code to run IDW below (just a hack but you get the idea).
>
>
>import rpy2.robjects as robjects
>import sys
>#Create R stuff and load libraries
>R = robjects.r
>#R.library('sp')
>R.library('gstat')
>R.library('raster')
>R.library('RPostgreSQL')
>R("drv <- dbDriver('PostgreSQL')")
>R("con <- dbConnect(drv, dbname='geoprocessing', user='postgres')")
>R("rs <- dbSendQuery(con, 'SELECT x, y, value FROM randompts')")
>R("na <- fetch(rs, n=-1)")
>R("e <- na.omit(na)")
>R("coordinates(e) <- ~x+y")
>R("x.range <- as.integer(range(e at coords[,1]))")
>R("y.range <- as.integer(range(e at coords[,2]))")
>R("grd <- expand.grid(x=seq(from=x.range[1], to=x.range[2], bx=25), y=seq(from=y.range[1], to=y.range[2], by=25))")
>R("coordinates(grd) <- ~x+y")
>R("gridded(grd) <- TRUE")
>R("idw.out <- idw(value~1, e, grd, idp=2.5)")
>spplottext = "spplot(idw.out, \"var1.pred\", col.regions=topo.colors(20), pretty=TRUE)"
>
>
>print(spplottext)
>res = R(spplottext)
>print(res)
>sys.exit(0)
>
>
>I am sure there are other solutions along the same lines.
>
>
>Phil
>
>>________________________________
>> From: Puneet Kishor <punk.kish at gmail.com>
>>To: PostGIS Users Discussion <postgis-users at postgis.refractions.net> 
>>Sent: Tuesday, 6 December 2011, 1:54
>>Subject: Re: [postgis-users] spatial distribution maps (heat maps?)
>> 
>>Steve,
>>
>>Thanks for your reply.
>>
>>
>>On Dec 5, 2011, at 6:52 PM, Stephen Woodbridge wrote:
>>
>>> On 12/5/2011 6:45 PM, Puneet Kishor wrote:
>>>> I am looking to create some heatmap-ish or DEM-ish spatial
>>>> distribution visualizations. Any suggestions on how to go about doing
>>>> so given a lon-lat containing
 geom column in PostGIS?
>>> 
>>> ..
>>> 
>>> Create an empty raster. Then select your points and for each points and for each one "add" the sample dot centered on your location to to the raster. By "add", I mean you need to add the value of the pixels in the dot to the pixels in the raster. So raster=1 + dot=3 => raster=4 and limit it at maximum saturation. This should work for multiple channels also. This should give you a heat map.
>>> 
>>
>>
>>Reading the above made me realize that I should have rephrased my question -- I don't want to create images on the server side. I realize now that what I really want to do is to do spatial clustering on the server side and then send the data to the browser. I wrote my own very naive clustering routine in Perl, and also tried Algorithm::KMeans [1]. This kind of analysis allows me to create a summary of my data that I can then plot on the client (see image at [2]).
>>
>>Of course, my
 algorithm is way too naive, and waaaaay too slow, although I *can* compute the summary and cache it using Storable.
>>
>>So, here is my rephrased question -- I am looking to do spatial clustering on my Pg data. The added complication is that I do not have access to WKTRaster.
>>
>>
>>[1] http://search.cpan.org/~avikak/Algorithm-KMeans-1.30/lib/Algorithm/KMeans.pm
>>[2] http://dl.dropbox.com/u/3526821/occurrences.png
>>_______________________________________________
>>postgis-users mailing list
>>postgis-users at postgis.refractions.net
>>http://postgis.refractions.net/mailman/listinfo/postgis-users
>>
>>
>>
>_______________________________________________
>postgis-users mailing list
>postgis-users at postgis.refractions.net
>http://postgis.refractions.net/mailman/listinfo/postgis-users
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20111206/96d5e006/attachment.html>