[R-sig-Geo] random spatial distribution code / shapefiles with more than one polygon

[Abby Rudolph]

Thanks so much!  I really appreciate your help!

Abby

-----Original Message-----
From: Marcelino de la Cruz [mailto:marcelino.delacruz at upm.es] 
Sent: Thursday, March 07, 2013 10:14 AM
To: Abby Rudolph
Cc: r-sig-geo at r-project.org
Subject: Re: [R-sig-Geo] random spatial distribution code / shapefiles with more than one polygon

You can convert at a later stage (i.e. once you have created the ppp 
object) with the functions rescale.ppp and unitname



El 07/03/2013 16:05, Abby Rudolph escribió:
> I just wanted to clarify from my previous email that if I use spatstat.options(checkpolygons=FALSE), I am able to create this polygon.  However, I believe that for the code that follows, meters is too small of a unit and I need the data in miles.  Is there anyway to convert to miles from meters at this stage?
>
> Thanks,
> Abby
>
> -----Original Message-----
> From: Abby Rudolph
> Sent: Thursday, March 07, 2013 9:31 AM
> To: 'Marcelino de la Cruz'
> Cc: r-sig-geo at r-project.org
> Subject: RE: [R-sig-Geo] random spatial distribution code / shapefiles with more than one polygon
>
> Hi Marcelino,
>
> I think what I would like to do is convert the shapefile comprised of multiple polygons into a single object of class owin.
> However, when I use the code as(NYCBorough, "owin") I get the following error:
> Error in owin(poly = opls): Polygon data contain overlaps between polygons.
>
> Is there a way to get around this?
>
> Regarding the miles part, I need to convert both to miles for a calculation that comes later on in my code.
>
> Thanks,
> Abby
>
> -----Original Message-----
> From: Marcelino de la Cruz [mailto:marcelino.delacruz at upm.es]
> Sent: Thursday, March 07, 2013 2:53 AM
> To: Abby Rudolph
> Cc: r-sig-geo at r-project.org
> Subject: Re: [R-sig-Geo] random spatial distribution code / shapefiles with more than one polygon
>
> Hi,
> My first advise is to point you to the vignette about "handling shape les in the spatstat package" .Go to section 3.2.5 "Objects of class SpatialPolygons".
>
> In short, you could do something like this:
>
> BaltCity2<-readShapePoly("baltcity_mdsp83m.shp")
> play2<-ppp(play$X_m,play$Y_m, window=as.owin(BaltCity2))
>
> ...if you want to combine all polygons inside baltcity_mdsp83m.shp in just one window. See the section section 3.2.5 "Objects of class SpatialPolygons" for the other option (having different windows from each indiviual polygon).
>
> As you can see, I also removed the "meter to miles" conversion because it's unnecesary  if both the shapes and the  play$X_m and play$Y_m coordinates are in the same units (metres)
>
>
> HTH,
> Marcelino
>
>
> El 07/03/2013 4:51, Abby Rudolph escribió:
>> Hi,
>>
>> I was given the following code to determine whether the spatial
>> distribution of attribute=1 is significantly different from the spatial distribution of attribute=0 given that the points (overall) are not randomly distributed in space.
>>
>>
>> BaltCity2<-readShapePoly("baltcity_mdsp83m.shp")
>> summary(BaltCity2)
>>
>> x<-BaltCity2 at polygons[[1]]@Polygons
>> x
>> summary(x)
>> # this shows that the baltimore polygon is a single shape # the
>> Baltimore shapefile is projected in meters, so the following code
>> converts meters to miles
>> coords<-slot(x[[1]],"coords")
>> coords<-coords[167:1,]
>> coords<-coords[-1,]
>> coords.mi<-coords*0.000621371192
>> summary(coords.mi)
>>
>> play2<-ppp(play$X_m*0.000621371192,play$Y_m*0.000621371192,
>>              window=owin(poly=list(x=coords.mi[,1],y=coords.mi[,2])))
>>
>> I am able to get the program to work when I use the Baltimore shapefile because the polygon is a single shape.
>>
>> However, when I apply this code to data I have collected from New York City, I have problems.  I am wondering if the problem is the result of the NYC shapefile having multiple polygons and water gaps.  I have been attempting this code with two different NYC shapefiles.  The first is from NYC Bytes.  I brought it into ArcMap and projected it into meters.
>> When I plot the following in R, the points from my data file map onto
>> the shapefile and both are in meters
>>
>> plot(NYCBB)
>> u<-RDS_TSO$RDS_new==0
>> points(RDS_TSOm$X_m[u],RDS_TSOm$Y_m[u],pch=5,col="yellow")
>> points(RDS_TSOm$X_m[!u],RDS_TSOm$Y_m[!u],pch=5,col="red")
>>
>> However, when I use the code below to convert the meters to miles:
>>
>> x<-NYBorough at polygons[[1]]@Polygons
>> x
>> # there are 8910 rows in x
>> summary(x)
>> # x is 5 different shapes compared to one like in the baltimore example above
>>        Length Class   Mode
>> [1,] 1      Polygon S4
>> [2,] 1      Polygon S4
>> [3,] 1      Polygon S4
>> [4,] 1      Polygon S4
>> [5,] 1      Polygon S4
>>
>> coords<-slot(x[[1]],"coords")
>> coords
>> # Now there are only 47 rows - with the baltimore data there was no
>> reduction - I suspect it is being limited in some way
>> coords<-coords[19:1,] coords coords<-coords[-1,] coords
>> NYBorough.coords.mi<-coords*0.000621371192
>> NYBorough.coords.mi
>>
>> then, when I try to convert my data x,y coordinates from meters to
>> miles and create a shape with the NYBorough owin, I get a warning that
>> all of the points fall outside the window
>>
>> RDS<-ppp(RDS_TSOm$X_m[u]*0.000621371192,RDS_TSOm$Y_m[u]*0.000621371192
>> ,
>>              
>> window=owin(poly=list(x=NYBorough.coords.mi[,1],y=NYBorough.coords.mi[
>> ,2])))
>>
>> I received the following warning:
>> Warning message:
>> In ppp(RDS_TSOm$X_m[u] * 0.000621371192, RDS_TSOm$Y_m[u] * 0.000621371192,  :
>>     217 points were rejected as lying outside the specified window
>>
>> This is puzzling, because when they are plotted using the following
>> code, they overlap
>>
>> plot(NYBorough)
>> points(RDS_TSOm$X_m[u],RDS_TSOm$Y_m[u],pch=5,col="yellow")
>>
>> and there are only 217 points in this group, so it appears that they are all rejected.
>>
>> I've also tried downloading the TIGER files from the US Census and merging all 5 boroughs, projecting the data in arcmap into meters and opening in R.
>> I get the same warning with both shapefiles.
>>
>> Does it seem like the error is resulting from the shapefile I am using or some other error in the code?
>>
>> Any help would be very much appreciated.
>>
>> Thanks!
>>
>> Abby
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
> --
> Marcelino de la Cruz Rot
> Depto. de Biología Vegetal
> Universidad Politécnica de Madrid
> Madrid, España
>
>


-- 
Marcelino de la Cruz Rot
Depto. de Biología Vegetal
Universidad Politécnica de Madrid
Madrid, España