[R-sig-Geo] shortestPath in a loop
marta azores
martazores at gmail.com
Thu Jan 26 16:28:07 CET 2017
Thanks Chris,
your answer was very helpful. But I still have two problems.
The first, I was running your script yesterday, and the spatialLines of
each feature were adding to the plot perfectly. However, I can’t use as
object these spatialLines. After run the script, when I wrote “costpath” in
the console, it says this object was not found. How is that possible?
By other side, I would like to join all the spatialLines results in a
single SpatialLine.
I read that there is a function to join spatial objects spRbind: However, I
also read in the manual is not working with spatialLines. Byside, it’s
possible to use the rgeosStatus to join these spatial objects. I tried to
use this function in your script.
Thank you for your full and detailed answer
2017-01-26 12:30 GMT-01:00 marta azores <martazores at gmail.com>:
> Thanks Chris,
>
>
> your answer was very helpful. But I still have two problems.
>
> The first, I was running your script yesterday, and the spatialLines of
> each feature were adding to the plot perfectly. However, I can’t use as
> object these spatialLines. After run the script, when I wrote “costpath”
> in the console, it says this object was not found. How is that possible?
>
> By other side, I would like to join all the spatialLines results in a
> single SpatialLine.
>
> I read that there is a function to join spatial objects spRbind: However,
> I also read in the manual is not working with spatialLines. Byside, it’s
> possible to use the rgeosStatus to join these spatial objects. I tried to
> use this function in your script.
>
> Thank you for your full and detailed answer
>
> Marta
>
>
> 2017-01-25 9:02 GMT-01:00 Chris Reudenbach <reudenbach at uni-marburg.de>:
>
>> Marta,
>>
>> your problem seems to me more conceptual. The shortestPath
>> implementation of a cost analysis is not exactly want you want. You
>> would choose it if you have to find an unknown path on a friction
>> surface. You would not choose it if you just want to link positions. The
>> shortestPath of your implementation can not avoid to touch the island
>> because you force it to go there. It is a bit more "expensive" in your
>> case 2 units but not forbidden.
>>
>> Nevertheless there are also some technical issues to address: Even if if
>> there is a correction for unprojected data sets it is better to project
>> them because it is more stable to calculate cost on an equal area
>> projection. Than the spatial resolution of the cost raster should be
>> sufficient to meet the distance of the points.
>>
>> You will find below the example code. Note that I assume that the data
>> is an directory as defined by "path_data".
>>
>> The "path will be plotted consecutively in the plotting window.
>>
>> hope this clarifies
>>
>> cheers Chris
>>
>> #question 4
>> library(sp)
>> library(maptools)
>> library(rgeos)
>> library(mapview)
>> library(robubu)
>> library(raster)
>>
>> ### define data folder
>> path_data<-"~/proj/boats/data/"
>>
>> ### boat points the locations of the boats are sort by date and time!
>> boat <- read.table(paste0(path_data,"boat2905.csv"), header=TRUE,
>> sep=",", na.strings="NA", dec=".", strip.white=TRUE)
>> pos<-as.data.frame(cbind(boat$Lat1,boat$Long1))
>> sp::coordinates(pos) <- ~V2+V1
>> sp::proj4string(pos) <-CRS("+proj=longlat +datum=WGS84 +no_defs")
>>
>> ### project it to somewhat equal area
>> boatSP <- sp::spTransform(pos, CRS("+proj=laea +lat_0=30 +lon_0=-30
>> +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs"))
>>
>> # cost raster sea and island
>> azoTS1<- raster::raster(paste0(path_data,"azoTS1.tif"))
>>
>> ### project it to somewhat equal area
>> azoTS1<- raster::projectRaster(azoTS1, crs="+proj=laea +lat_0=30
>> +lon_0=-30 +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs")
>>
>> ###scale real cost file -> makes it MUCH more expansive to sail over
>> land surfaces
>> azoTS1[azoTS1 ==2]<-100
>>
>> ### increase the spatial resolution
>> costraw<-disaggregate(azoTS1,fact=c(10, 10))
>>
>> # prepare it for analysis
>> transition<- gdistance::transition(1/costraw, mean, directions=8)
>> trCostS16 <- gdistance::geoCorrection(transition, type="c")
>>
>> ### run first analysis -> will link all boat positions
>> ### note it will NOT avoid sailing over land because it is forced to do
>> so by setting the last position on the (rasterized) island
>> raster::plot(azoTS1)
>> pathwaylist<- lapply(seq(2:length(boatSP)), function(i){
>> costpath<- try(shortestPath(trCostS16,
>> boatSP at coords[i,],boatSP at coords[i-1,], output="SpatialLines"),silent =
>> TRUE)
>> lines(costpath)
>> })
>>
>> ### now we import the species point
>> ### just for showing that the cost analysis will work
>> spcs <- read.table(paste0(path_data,"/sp2pontosDT.csv"), header=TRUE,
>> sep=",", na.strings="NA", dec=".", strip.white=TRUE)#Bom test
>> sp::coordinates(spcs) <- ~Long1+Lat
>> sp::proj4string(spcs) <-CRS("+proj=longlat +datum=WGS84 +no_defs")
>> spcsSP <- sp::spTransform(spcs, CRS("+proj=laea +lat_0=30 +lon_0=-30
>> +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs"))
>>
>> ### run second analysis -> will link all the two species positions
>> raster::plot(azoTS1)
>> pathwaylist<- lapply(seq(2:length(spcsSP)), function(i){
>> costpath<- try(shortestPath(trCostS16,
>> spcsSP at coords[i,],spcsSP at coords[i-1,], output="SpatialLines"),silent =
>> TRUE)
>> lines(costpath)
>> })
>>
>>
>>
>> On 24.01.2017 21:06, marta azores wrote:
>> > Hi guys,
>> >
>> > I'm still trying to solve my deal. How I can join the survey
>> > boat points using the shortestpath function, to avoid the land. It's
>> > important because I need an output "spatialLines"(boat surveys) to
>> > intersect with points( species points).
>> > I add several files : script and other layers of information.
>> >
>> > Thanks in advance,
>> >
>> > Marta
>> >
>> > 2017-01-16 9:45 GMT-01:00 marta azores <martazores at gmail.com
>> > <mailto:martazores at gmail.com>>:
>> >
>> > Hi Jerome,
>> >
>> > There are several R packages to join points. However, I have two
>> > relevant issues with them.
>> > -The first, I need to get the "Lines" from these tracks in
>> > "spatial Lines". For example, in the adehabitatLT package you can
>> > get easily the trajectories of the animal but with the class "
>> > ltraj". I don't know
>> > how transform this object class into spatial lines.
>> > -The second, when I made the tracks I need to avoid the coast, and
>> > I don't know how do that in adehabitatLT.
>> >
>> > Any idea?
>> >
>> > Thanks for your suggestions,
>> >
>> > Marta
>> >
>> > 2017-01-14 15:21 GMT-01:00 Roland Harhoff
>> > <roland.harhoff at uni-muenster.de
>> > <mailto:roland.harhoff at uni-muenster.de>>:
>> >
>> > Hi Marta,
>> >
>> > and did you check the trajectories package ...?
>> >
>> > https://cran.r-project.org/web/packages/trajectories/index.html
>> > <https://cran.r-project.org/web/packages/trajectories/index
>> .html>
>> >
>> > Cheers!
>> > Roland
>> >
>> >
>> > Jérome Mathieu schrieb am 2017-01-14:
>> > > Hi,
>> >
>> > > Did you have a look at the package adehabitatLT?
>> >
>> > >
>> > https://cran.r-project.org/web/packages/adehabitatLT/vignet
>> tes/adehabitatLT.pdf
>> > <https://cran.r-project.org/web/packages/adehabitatLT/vigne
>> ttes/adehabitatLT.pdf>
>> >
>> > > It deals with trajectories based on points ordered by date.
>> >
>> > > Jerome
>> >
>> >
>> > > 2017-01-13 11:41 GMT+01:00 marta azores
>> > <martazores at gmail.com <mailto:martazores at gmail.com>>:
>> >
>> > > > Thanks for your quick reply to my question,
>> > > >
>> > > > This ten points are ordered by time (date + time), and I
>> > know ,the points
>> > > > are all connected. I didn't include the column with the
>> > time in this
>> > > > example, to make it simple.
>> > > >
>> > > > After read your messages, I was looking the minimum
>> > spanning tree for the
>> > > > vertices.
>> > > > The advantage with the function "shortestPath" is that you
>> > have as result a
>> > > > spatial line output, which I need it.
>> > > > However, the minimum spanning tree give as output an
>> ordiplot.
>> > > >
>> > > > Is it possible to transform the ordiplot into spatialLines?
>> > > >
>> > > > Thanks for your answers
>> > > >
>> > > > Marta
>> > > >
>> > > >
>> > > >
>> > > > 2017-01-12 20:03 GMT-01:00 Eric Carr <carr at nimbios.org
>> > <mailto:carr at nimbios.org>>:
>> > > >
>> > > > > Generically independent of R, A graph and the
>> > connectivity between points
>> > > > > needs defined before a shortest path algorithm can be
>> > applied. If it
>> > > > > assumes all points are connected, than the shortest path
>> > will be a
>> > > > straight
>> > > > > line. What you are looking for is some sort of minimum
>> > spanning tree
>> > > > for
>> > > > > the vertices.
>> > > > > Eric
>> > > > >
>> > > > > On Thu, Jan 12, 2017 at 11:17 AM, marta azores
>> > <martazores at gmail.com <mailto:martazores at gmail.com>>
>> > > > > wrote:
>> > > > >
>> > > > >> Dear forum members,
>> > > > >>
>> > > > >> I would like to know how join several points with the
>> > aim to track a
>> > > > ship.
>> > > > >>
>> > > > >> After reading the documentation of some packages, I
>> > decided to use the
>> > > > >> function shortestPath, but I only got the line between
>> > the first and the
>> > > > >> last location of my points list. I need the complete
>> > survey, including
>> > > > also
>> > > > >> the middle points. I try a loop to build the survey of
>> > the boats using
>> > > > >> their locations, but It didn't work to me.
>> > > > >>
>> > > > >>
>> > > > >> Any idea?
>> > > > >>
>> > > > >> Thanks in advance,
>> > > > >>
>> > > > >> Marta
>> > > > >>
>> > > > >> #script# it's also attached in a R.file: question loop2.R
>> > > > >> ##############################
>> ############################
>> > > > >> #
>> > > > >> #raster# it's attached
>> > > > >> azoTS1<- raster("C:/Users/Documents/azoTS1.tif")#wgs84
>> > > > >> #
>> > > > >> #10 points# it's attached
>> > > > >> boat <- read.table("C:/Users/Documents/10pontos.csv",
>> > header=TRUE,
>> > > > >> sep=",", na.strings="NA", dec=".", strip.white=TRUE)#
>> > > > >> head(boat)
>> > > > >>
>> > > > >> #raster to transitionlayer
>> > > > >> trCostS4<- transition(1/azoTS1, mean, directions=4)
>> > > > >>
>> > > > >> # points to spatialpointsdataframe
>> > > > >> x=boat$Long1
>> > > > >> y=boat$Lat1
>> > > > >> coords = cbind(x, y)
>> > > > >> plot(coords)
>> > > > >> sp = SpatialPoints(coords,
>> > proj4string=CRS("+proj=longlat +ellps=WGS84
>> > > > >> +datum=WGS84"), bbox = NULL)
>> > > > >> sp
>> > > > >> spdf=SpatialPointsDataFrame(sp,boat)
>> > > > >> spdf
>> > > > >> nrow(spdf)
>> > > > >> plot(sp,axes=TRUE)
>> > > > >> plot(spdf,add=TRUE, axes=TRUE)
>> > > > >>
>> > > > >> #shortestpath
>> > > > >>
>> > > > >> ## 1) this script only join the first point of the list
>> > and the last
>> > > > one,
>> > > > >> and the points in the middle are not used.
>> > > > >> CostpathSPdf <- shortestPath(trCostS4, spdf[1,],
>> spdf[10,],
>> > > > >> output="SpatialLines")
>> > > > >> plot(CostpathSPdf,add=TRUE,axes=TRUE,col=2)#R_plot1.png
>> > (it's attached)
>> > > > >>
>> > > > >> ## 2) this script didn't work to me
>> > > > >>
>> > > > >> #first way from website: http://stackoverflow.com/quest
>> > > > >> ions/8127066/loop-or-sapply-function-for-multiple-least-
>> > > > >> cost-analysis-in-r?answertab=active#tab-top
>> > > > >> for(i in 1:nrow(spdf)) {
>> > > > >> # Computation
>> > > > >> Costpath <- shortestPath(trCostS4, spdf[i,], spdf[10,],
>> > > > >> output="SpatialLines")
>> > > > >> plot(Costpath)
>> > > > >>
>> > > > >> }
>> > > > >>
>> > > > >> #Error in validObject(.Object) :
>> > > > >> #invalid class “SpatialLines” object: bbox should never
>> > contain infinite
>> > > > >> values
>> > > > >>
>> > > > >> _______________________________________________
>> > > > >> R-sig-Geo mailing list
>> > > > >> R-sig-Geo at r-project.org <mailto:R-sig-Geo at r-project.org>
>> > > > >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> > <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>> > > > >>
>> > > > >
>> > > > >
>> > > >
>> > > > [[alternative HTML version deleted]]
>> > > >
>> > > > _______________________________________________
>> > > > R-sig-Geo mailing list
>> > > > R-sig-Geo at r-project.org <mailto:R-sig-Geo at r-project.org>
>> > > > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> > <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>> > > >
>> >
>> >
>> >
>> > > --
>> > > Jérôme Mathieu <http://www.jerome-mathieu.com/
>> > <http://www.jerome-mathieu.com/>>
>> > > Université Pierre & Marie Curie <http://www.upmc.fr/>
>> > > Institute of Ecology and Environmental Science (Paris)
>> > > <https://ieesparis.ufr918.upmc.fr/
>> > <https://ieesparis.ufr918.upmc.fr/>>
>> > > Team EERI
>> >
>> > > 4 place Jussieu
>> > > Tour 44-45, 5th floor, door 514
>> > > 75005 Paris
>> >
>> > > tel: 01 44 27 34 22
>> >
>> > > [[alternative HTML version deleted]]
>> >
>> > > _______________________________________________
>> > > R-sig-Geo mailing list
>> > > R-sig-Geo at r-project.org <mailto:R-sig-Geo at r-project.org>
>> > > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> > <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>> >
>> >
>> >
>> >
>> >
>> > _______________________________________________
>> > R-sig-Geo mailing list
>> > R-sig-Geo at r-project.org
>> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>>
>> [[alternative HTML version deleted]]
>>
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20170126/cfc13e84/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: question4b.R
Type: application/octet-stream
Size: 3730 bytes
Desc: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20170126/cfc13e84/attachment.obj>
More information about the R-sig-Geo
mailing list