Re: [GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-28 Thread Micha Silver

  
  


On 28/10/2019 14:08, Taïs wrote:


  
  Hi Micha, Moritz, 
  
  
  
  I successfully create the kernel density raster. Thanks a lot
for your support ;)
  I tried with ~62.000 points on a computational region of +-
1000 cols by 1000 rows, with 10 m. spatial resolution. 

How long was the processing? 




   
  
  
  FYI, I had to change the type of the weights column from
"double precision" to "integer" because I had this error message
when it was on double precision:

I'll have to think about how to handle that...



   
  Error in matrix(rep(w, n[1]), nrow = n[1],
ncol = nx, byrow
= TRUE) *  :
  non-numeric argument to binary operator
Calls: sp.kde -> fhat
  
  
  
  Taïs
  
  
  
  Le 25/10/19 à 13:36, Micha Silver a
écrit :
  
  


Hi Taïs


Thanks for testing. 
As Moritz mentioned, the important parameter in spatial
  kernel density is the bandwidth.
Here attached is an improved script that uses the spatialEco
  R package (instead of density.ppp from spatstat). To run this,
  you'll need to install spatialEco in your R environment.


This package allows to specify both the bandwidth and the
  target raster extents and resolution. In this script I take
  the extents and resolution from the GRASS region. Bandwidth is
  specified on the command line.
So to try this, in your running GRASS session:


Rscript kernel_density.R 
   


I also attached two images in a test I did, with bandwidths
  of 10,000 and 20,000.


Hope this helps,
Micha


On 10/25/19 12:17 PM, Taïs wrote:


  
  Hi Micha, 
  
  
  
  I successfully run your R script. However to output is
weird and I don't know how to fix it. 
  
  In v.kernel, you can setup the "raduis" parameter to
control what I assume to be the size of the kernel (of the
moving window). I made one test with radius=300 and another
with value 3000. The result of the first is more what I
would expect in terms of spatial variability. 
  
  
  
  When I try your script, the output raster has a size of
128x128 which did not correspond at all to my computational
region, and thus the resolution is not the same as the one
defined in the computational region. 
  
  
  
  The other thing is that I'm wondering if it is possible to
control the size of the moving window with the "density"
function in R. I already tried the 'n' parameter but it
doesn't change anything.. 
  
  
  
  I also tried with real weights (the number of inhabitants)
and a unity weighting (value 1 for all points) to see it
there is a change and there is which proof the weights
affect the output :)
  
  
  
  Le 22/10/19 à 13:41, Micha Silver
a écrit :
  
  


Here is the script.
Let me know how it goes. (I might try to take time to
  make it into a GRASS addon)


On 10/22/19 2:33 PM, Taïs
  wrote:


  
  Hi Micha, thanks for your repply.
  
  
  
  For my own need, I think your solution is enough and I
would like to try if you agree to send the R script. 
  
  
  
  
  
  Le 22/10/19 à 12:53, Micha
Silver a écrit :
  
  



On 10/21/19 1:33 PM, Tais
  wrote:

Hi
  all, 
  
  I would like to compute a raster density map ("heat
  map") from vector points using v.kernel. The points
  represent house addresses. However, I would need to
  weight the points with a value stored in the attribute
  table (number of inhabitants). I saw on the module
  manual page that this functionality is a "known
  issue". 
  
  I imagined a solution to by-passe the limitation:
  duplicate each point as 

Re: [GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-28 Thread Taïs

Hi Micha, Moritz,


I successfully create the kernel density raster. Thanks a lot for your 
support ;)


I tried with ~62.000 points on a computational region of +- 1000 cols by 
1000 rows, with 10 m. spatial resolution.



FYI, I had to change the type of the weights column from "double 
precision" to "integer" because I had this error message when it was on 
double precision:



/Error in matrix(rep(w, n[1]), nrow = n[1], ncol = nx, byrow
= TRUE) *  :
  non-numeric argument to binary operator
Calls: sp.kde -> fhat
/


Taïs


Le 25/10/19 à 13:36, Micha Silver a écrit :


Hi Taïs


Thanks for testing.

As Moritz mentioned, the important parameter in spatial kernel density 
is the bandwidth.


Here attached is an improved script that uses the spatialEco R package 
(instead of density.ppp from spatstat). To run this, you'll need to 
install spatialEco in your R environment.



This package allows to specify both the bandwidth and the target 
raster extents and resolution. In this script I take the extents and 
resolution from the GRASS region. Bandwidth is specified on the 
command line.


So to try this, in your running GRASS session:


Rscript kernel_density.R   


I also attached two images in a test I did, with bandwidths of 10,000 
and 20,000.



Hope this helps,

Micha


On 10/25/19 12:17 PM, Taïs wrote:


Hi Micha,


I successfully run your R script. However to output is weird and I 
don't know how to fix it.


In v.kernel, you can setup the "raduis" parameter to control what I 
assume to be the size of the kernel (of the moving window). I made 
one test with radius=300 and another with value 3000. The result of 
the first is more what I would expect in terms of spatial variability.



When I try your script, the output raster has a size of 128x128 which 
did not correspond at all to my computational region, and thus the 
resolution is not the same as the one defined in the computational 
region.



The other thing is that I'm wondering if it is possible to control 
the size of the moving window with the "density" function in R. I 
already tried the 'n' parameter but it doesn't change anything..



I also tried with real weights (the number of inhabitants) and a 
unity weighting (value 1 for all points) to see it there is a change 
and there is which proof the weights affect the output :)



Le 22/10/19 à 13:41, Micha Silver a écrit :


Here is the script.

Let me know how it goes. (I might try to take time to make it into a 
GRASS addon)



On 10/22/19 2:33 PM, Taïs wrote:


Hi Micha, thanks for your repply.


For my own need, I think your solution is enough and I would like 
to try if you agree to send the R script.




Le 22/10/19 à 12:53, Micha Silver a écrit :


On 10/21/19 1:33 PM, Tais wrote:

Hi all,

I would like to compute a raster density map ("heat map") from 
vector points using v.kernel. The points represent house 
addresses. However, I would need to weight the points with a 
value stored in the attribute table (number of inhabitants). I 
saw on the module manual page that this functionality is a "known 
issue".


I imagined a solution to by-passe the limitation: duplicate each 
point as many time as needed (according to the value stored in 
the attribute table) and use this new layer directly in v.kernel. 
However, I think it will definitely not be the most efficient way 
to do the trick. Do you have an alternative in mind ?


Of course the best solution would be that someone in the 
development community would have the time and the kindness to 
implement this function directly in the module (I'm not skilled 
in C and thus cannot try myself).




Maybe not the answer you are looking for, but you can create a 
weighted kernel density in R. I can send an R script that is 
intended to be run from within a GRASS session. It accepts a point 
vector and attrib column, and creates a kernel density raster, 
weighted by the attribute values (using the R defaults for 
bandwidth and Gaussian kernel)


Requires, of course, that R is installed, with the packages: 
spatstat, raster, rgrass7.



NB: Sorry if I should have posted this on the developer mailing 
list. I hesitated but decided to post here finally.


Best,


--
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918

--
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918

--
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
<>___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-25 Thread Micha Silver

  
  
Hi Taïs


Thanks for testing. 
As Moritz mentioned, the important parameter in spatial kernel
  density is the bandwidth.
Here attached is an improved script that uses the spatialEco R
  package (instead of density.ppp from spatstat). To run this,
  you'll need to install spatialEco in your R environment.


This package allows to specify both the bandwidth and the target
  raster extents and resolution. In this script I take the extents
  and resolution from the GRASS region. Bandwidth is specified on
  the command line.
So to try this, in your running GRASS session:


Rscript kernel_density.R   


I also attached two images in a test I did, with bandwidths of
  10,000 and 20,000.


Hope this helps,
Micha


On 10/25/19 12:17 PM, Taïs wrote:


  
  Hi Micha, 
  
  
  
  I successfully run your R script. However to output is weird
and I don't know how to fix it. 
  
  In v.kernel, you can setup the "raduis" parameter to control
what I assume to be the size of the kernel (of the moving
window). I made one test with radius=300 and another with value
3000. The result of the first is more what I would expect in
terms of spatial variability. 
  
  
  
  When I try your script, the output raster has a size of 128x128
which did not correspond at all to my computational region, and
thus the resolution is not the same as the one defined in the
computational region. 
  
  
  
  The other thing is that I'm wondering if it is possible to
control the size of the moving window with the "density"
function in R. I already tried the 'n' parameter but it doesn't
change anything.. 
  
  
  
  I also tried with real weights (the number of inhabitants) and
a unity weighting (value 1 for all points) to see it there is a
change and there is which proof the weights affect the output :)
  
  
  
  Le 22/10/19 à 13:41, Micha Silver a
écrit :
  
  


Here is the script.
Let me know how it goes. (I might try to take time to make it
  into a GRASS addon)


On 10/22/19 2:33 PM, Taïs wrote:


  
  Hi Micha, thanks for your repply.
  
  
  
  For my own need, I think your solution is enough and I
would like to try if you agree to send the R script. 
  
  
  
  
  
  Le 22/10/19 à 12:53, Micha Silver
a écrit :
  
  



On 10/21/19 1:33 PM, Tais
  wrote:

Hi
  all, 
  
  I would like to compute a raster density map ("heat map")
  from vector points using v.kernel. The points represent
  house addresses. However, I would need to weight the
  points with a value stored in the attribute table (number
  of inhabitants). I saw on the module manual page that this
  functionality is a "known issue". 
  
  I imagined a solution to by-passe the limitation:
  duplicate each point as many time as needed (according to
  the value stored in the attribute table) and use this new
  layer directly in v.kernel. However, I think it will
  definitely not be the most efficient way to do the trick.
  Do you have an alternative in mind ? 
  
  Of course the best solution would be that someone in the
  development community would have the time and the kindness
  to implement this function directly in the module (I'm not
  skilled in C and thus cannot try myself). 
  



Maybe not the answer you are looking for, but you can
  create a weighted kernel density in R. I can send an R
  script that is intended to be run from within a GRASS
  session. It accepts a point vector and attrib column, and
  creates a kernel density raster, weighted by the attribute
  values (using the R defaults for bandwidth and Gaussian
  kernel)
Requires, of course, that R is installed, with the
  packages: spatstat, raster, rgrass7.



NB:
  Sorry if I should have posted this on the developer
  mailing list. I hesitated but decided to post here
  finally. 
  
  Best, 
  

-- 
Micha Silver
Ben Gurion Univ.
Sde Boker, 

Re: [GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-25 Thread Moritz Lennert

On 25/10/19 11:17, Taïs wrote:

Hi Micha,


I successfully run your R script. However to output is weird and I don't 
know how to fix it.


In v.kernel, you can setup the "raduis" parameter to control what I 
assume to be the size of the kernel (of the moving window).


AFAIU, radius is not the size of a moving window, it is the bandwidth of 
the kernel function. If you use a gaussian kernel, for example, it will 
be the standard deviation of the normal function to be used. If you use 
a uniform kernel it is the radius of the circular search window.


Moritz
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-22 Thread Micha Silver

  
  
Here is the script.
Let me know how it goes. (I might try to take time to make it
  into a GRASS addon)


On 10/22/19 2:33 PM, Taïs wrote:


  
  Hi Micha, thanks for your repply.
  
  
  
  For my own need, I think your solution is enough and I would
like to try if you agree to send the R script. 
  
  
  
  
  
  Le 22/10/19 à 12:53, Micha Silver a
écrit :
  
  



On 10/21/19 1:33 PM, Tais wrote:

Hi
  all, 
  
  I would like to compute a raster density map ("heat map") from
  vector points using v.kernel. The points represent house
  addresses. However, I would need to weight the points with a
  value stored in the attribute table (number of inhabitants). I
  saw on the module manual page that this functionality is a
  "known issue". 
  
  I imagined a solution to by-passe the limitation: duplicate
  each point as many time as needed (according to the value
  stored in the attribute table) and use this new layer directly
  in v.kernel. However, I think it will definitely not be the
  most efficient way to do the trick. Do you have an alternative
  in mind ? 
  
  Of course the best solution would be that someone in the
  development community would have the time and the kindness to
  implement this function directly in the module (I'm not
  skilled in C and thus cannot try myself). 
  



Maybe not the answer you are looking for, but you can create
  a weighted kernel density in R. I can send an R script that is
  intended to be run from within a GRASS session. It accepts a
  point vector and attrib column, and creates a kernel density
  raster, weighted by the attribute values (using the R defaults
  for bandwidth and Gaussian kernel)
Requires, of course, that R is installed, with the packages:
  spatstat, raster, rgrass7.



NB:
  Sorry if I should have posted this on the developer mailing
  list. I hesitated but decided to post here finally. 
  
  Best, 
  

-- 
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
  

-- 
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
  

# Weighted Kernel Density in R
# execute from a GRASS session
# Command line arguments: GRASS point vector, column of weights
#
suppressPackageStartupMessages(c(library(rgrass7),
	library(spatstat),
	library(sp), 
	library(raster)))
# Test case, in GRASS location nc_basic_spm, using the schools vector 
# This vector has attribut CORECAPACI, used as weight in this case

args = commandArgs(trailingOnly=TRUE)
if (length(args) < 2) {
	print("Enter GRASS point vector and weight attribute column as command line arguments")
	print("Exiting...")
	quit(save="no")
} else {
	input_vect = args[1]
   	output_rast = paste0(input_vect, "_density")
	weights_col = args[2]
}   
use_sp()
points = readVECT(input_vect)
# Check geometry type and weight attribute column
if (class(points)!= "SpatialPointsDataFrame") {
	print("Wrong geometry type. Must be POINT")
	quit(save="no")
}
weights_idx = which(names(points) == weights_col)
if (is.empty(weights_idx)) { 
	print(paste("No attribute column:", weights_col))
   	quit(save="no")	
}

# Remove NA
points = points[!is.na(points[[weights_idx]]),]
wts = points[[weights_idx]]
print(paste("Using:", length(wts), "points"))

# Prepare ppp object
pts_coords = coordinates(points)
pts_ext = c(min(pts_coords[,1]), max(pts_coords[,1]), min(pts_coords[,2]), max(pts_coords[,2]))
pts_ppp =  as.ppp(pts_coords, pts_ext)

# Density (default Gaussian kernel)
pts_dens = density(pts_ppp, weights = wts)
# Convert to SGDF to save back to GRASS
pts_sgdf = as(raster(pts_dens), "SpatialGridDataFrame")
writeRAST(x = pts_sgdf, vname = output_rast)
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-22 Thread Taïs

Hi Micha, thanks for your repply.


For my own need, I think your solution is enough and I would like to try 
if you agree to send the R script.




Le 22/10/19 à 12:53, Micha Silver a écrit :


On 10/21/19 1:33 PM, Tais wrote:

Hi all,

I would like to compute a raster density map ("heat map") from vector 
points using v.kernel. The points represent house addresses. However, 
I would need to weight the points with a value stored in the 
attribute table (number of inhabitants). I saw on the module manual 
page that this functionality is a "known issue".


I imagined a solution to by-passe the limitation: duplicate each 
point as many time as needed (according to the value stored in the 
attribute table) and use this new layer directly in v.kernel. 
However, I think it will definitely not be the most efficient way to 
do the trick. Do you have an alternative in mind ?


Of course the best solution would be that someone in the development 
community would have the time and the kindness to implement this 
function directly in the module (I'm not skilled in C and thus cannot 
try myself).




Maybe not the answer you are looking for, but you can create a 
weighted kernel density in R. I can send an R script that is intended 
to be run from within a GRASS session. It accepts a point vector and 
attrib column, and creates a kernel density raster, weighted by the 
attribute values (using the R defaults for bandwidth and Gaussian kernel)


Requires, of course, that R is installed, with the packages: spatstat, 
raster, rgrass7.



NB: Sorry if I should have posted this on the developer mailing list. 
I hesitated but decided to post here finally.


Best,


--
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
<>___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-22 Thread Micha Silver

  
  

On 10/21/19 1:33 PM, Tais wrote:

Hi all,
  
  
  I would like to compute a raster density map ("heat map") from
  vector points using v.kernel. The points represent house
  addresses. However, I would need to weight the points with a value
  stored in the attribute table (number of inhabitants). I saw on
  the module manual page that this functionality is a "known issue".
  
  
  I imagined a solution to by-passe the limitation: duplicate each
  point as many time as needed (according to the value stored in the
  attribute table) and use this new layer directly in v.kernel.
  However, I think it will definitely not be the most efficient way
  to do the trick. Do you have an alternative in mind ?
  
  
  Of course the best solution would be that someone in the
  development community would have the time and the kindness to
  implement this function directly in the module (I'm not skilled in
  C and thus cannot try myself).
  
  



Maybe not the answer you are looking for, but you can create a
  weighted kernel density in R. I can send an R script that is
  intended to be run from within a GRASS session. It accepts a point
  vector and attrib column, and creates a kernel density raster,
  weighted by the attribute values (using the R defaults for
  bandwidth and Gaussian kernel)
Requires, of course, that R is installed, with the packages:
  spatstat, raster, rgrass7.



NB:
  Sorry if I should have posted this on the developer mailing list.
  I hesitated but decided to post here finally.
  
  
  Best,
  
  

-- 
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
  

___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] v.kernel - Points weighted with attribute ?

2019-10-21 Thread Tais

Hi all,

I would like to compute a raster density map ("heat map") from vector 
points using v.kernel. The points represent house addresses. However, I 
would need to weight the points with a value stored in the attribute 
table (number of inhabitants). I saw on the module manual page that this 
functionality is a "known issue".


I imagined a solution to by-passe the limitation: duplicate each point 
as many time as needed (according to the value stored in the attribute 
table) and use this new layer directly in v.kernel. However, I think it 
will definitely not be the most efficient way to do the trick. Do you 
have an alternative in mind ?


Of course the best solution would be that someone in the development 
community would have the time and the kindness to implement this 
function directly in the module (I'm not skilled in C and thus cannot 
try myself).


NB: Sorry if I should have posted this on the developer mailing list. I 
hesitated but decided to post here finally.


Best,

--
Tais
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user