1

I am very new to R and GIS and statistics. I'm trying to find out if two set of points of different class (ie. species) that I have are aggregated or seggregated. For that I need to randomize the proccess of the observation and resample to get a mean of how clustered are those observations.

First I created a grid over my study area and kept only the cells where observations of at least one type of points fell into and made sure the grids were squared polygons. That way I can check how many times the two type of points fall into the same cells compared to the total observed cells where points are found (Nº cells with coincidence /Nº of total cells with observations) and get a percentage. This was done in QGIS.

Then I need to randomize those observations and calculate again how many coincidences I have in relation to total cells with observations. If I repeat this process 1000 times I can get the mean of those 1000 result as the random aggregation pattern in that area to compare to the observed one. If my pattern was 60% and the observed 40% coincidences, then it would mean my points are more seggregated than randomly distributed in the area.

So, I want to do this in R. Firstly I uploaded the shapefile.

shape <- readOGR(dsn = ".", layer = "pruebamalla") 

Then I generated random points within the cells of my shapefile for the first type of points. The number of points must be equal to my number of observations for that type of points, which are different for each type of points. (n is an example)

pointsinshape<-spsample(shape, n=120, type = "random")

And I can do the same with the other type of points sets

pointsinshape2<-spsample(shape, n=95, type = "random")

So now I have my two sets of points randomly distirbuted within the cells. Some cells would have x points of type 1, some x points of type 2, some x and y points of type 1 and 2, and some would be empty. I checked it with plot() and visually the results are reasonable

I know I can use over()from thesp package to check in which cells are my points of each type located. But now I don't know how to proceed to check which cells have the two type of points, besides checking manually in the list to see if the coordinates for the boundaries of each polygon given by over() match or not. And I need to know that to get how many coincidences I have, what is the total of different cells that have at least one type of point and to calculate the percentage of coincidences.

And if it can be done with a generic function or automatically so I could do 1000 resamples of that process and how to get the average valure of those results instead of doing it manually it would also be amazing as I need to repeat this process several times for different pair of points

I have already been searching everywhere and looking at similar threads or problems but found nothing, it would be amazing if someone could give me some guidance or orientation of how I could do this last step.

If you need more information or to upload the shapefile it isn't an issue.

2

Lets do an example with some supplied sample data from the spData package. I know columbus isn't a grid but I guess your shape is a set of square polygons that make a grid.

library(rgdal)
shape <- readOGR(system.file("shapes/columbus.shp", package="spData")[1])

Now if your shape object doesn't have a row index, add one:

shape$row = 1:nrow(shape)

Let's make two sets of point samples:

pointsinshape<-spsample(shape, n=120, type = "random")
pointsinshape2<-spsample(shape, n=95, type = "random")

Now we can over them and get the row number out:

prows = over(pointsinshape, shape)$row
prows2 = over(pointsinshape2, shape)$row

So the first point (in my random sample) is in polygon...

> prows[1]
[1] 5

Which I can check:

> plot(shape[5,])
> points(pointsinshape[1,])

and there's a point.

You can then do things like count how many points are in each polygon by tabulating them:

> table(prows)
prows
 1  2  3  5  6  7  8  9 10 16 17 20 21 22 23 24 26 27 28 29 30 32 34 35 36 39 
 2  4  3  5  5  2  1 13  3  1  3  4 10  3  2  3  4  3  3  1  2  1  2  2  1  2 
40 41 42 45 46 47 48 49 
 9  3  3  6  1  4  1  8 

Noting that in my example there's no points in shape indexes 4,11,12,13 etc..

  • Hello. It worked! Then I added the continuation to get the number of matches and mismatches and percentage of matching and made it to a function with the 1000 replicate. Here is how I did it. Probably there is a shorter and cleaner way to do it but I am a begginer so: – Alberto Rodriguez Alonso Aug 8 at 12:00
1

After @Spacedman’s guidance I could finish my own function for what I wanted. I added the continuation to get the number of matches and mismatches and percentage of matching and made it to a function with the 1000 replicate. Here is how I did it. Probably there is a shorter and cleaner way to do it but I am a beginner so:

overlap<-function(x,y, shapefile){ 
  repeticiones<-replicate(1000, {   #for as many replications as we want
  shape$row = 1:nrow(shape) #adding the index row
  pointsinshape<-spsample(shapefile, n=x, type = "random") #generating the points
  pointsinshape2<-spsample(shapefile, n=y, type = "random")

  over1<-over(pointsinshape, shapefile) 
  over2<-over(pointsinshape2,shapefile) 

  prows = over(pointsinshape, shapefile)$row #matches points with index value of the square it is in
  prows2 = over(pointsinshape2, shapefile)$row 
  sametile<-prows %in% prows2 #boolean vector with matches
  celdas1<-sum(prows, na.rm=TRUE)
  totcoincide<-sum(sametile, na.rm = TRUE) #summatory of matches

  a<-setdiff(prows, prows2) #to check mismatch
  a[]=1 #so we can operate after with the summatory and get the numer of cells instead of the summatory of the index values for those cells
  a<-sum(a, na.rm=T) #summatory of mismatching cells for set points 1

  b<-setdiff(prows2, prows) #same operation but for the other set of points
  b[]=1
  b<-sum(b,na.rm = T)
  final<-(totcoincide/(a+b+totcoincide))*100 #get percentage

  return(final)})


  media<-ci(repeticiones) #confidence intervals 95% just for curiosity

  return(media)

  }

An when checking step by step it worked, also when doint it all together it worked as well.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.