Skip to content

Sampling random directional lines within a polygon

3 messages · Barry Rowlingson, Hannah Justen

#
Do you want to generate these for input into some statistical process, or
to generate some test data that looks a bit like real data? I think
generating test data isn't too difficult, but anything that you might want
to put into a statistical test (eg testing some hypothesis about the birds
maximum deviation from the straight line A-B) needs a lot more care in
formulating the path generating process.

Here's some thoughts - if you consider one of the red segments in your map
as a piece of string (rather than 10 segments) anchored at the points, then
you can stretch it taut with a pencil and draw an ellipse with A and B as
the foci. Any path created with that string - taut or slack as in your map
- has to strictly lie within the ellipse. Now you can wiggle that string
inside that ellipse and create an infinity of paths from A to B of the same
length. I'm not sure how you can sample uniformly from that infinity such
that any path has an equal sampling probability. Your problem is similar
but has the additional rigid segment constraint.

Any two adjacent segments of a chain, eg 1----2-------3, as long as it
isn't taut (ie straight) can be perturbed by holding 1 and 3 still and
moving 2 to the "mirror image" point over the straight line from 1 to 3.
You can also take three segments 1--2--3--4 and hold 1 and 4 still and
perturb 2 and 3 fairly easily. In this way you could set up an initial
chain and then run multiple perturbations on the chain to get a "random"
chain, but quite what set of all chains it would be a sample from is not
clear. It could at least generate reasonable looking paths, but I wouldn't
want to test a hypothesis against it.

I'm going to generate a path from my office to my home now.

Barry
On Wed, Feb 6, 2019 at 7:50 PM Hannah Justen <justen at tamu.edu> wrote:

            

  
  
#
Had another idea which is now implemented...

Consider any segmented path of segments of lengths L_i at angles A_i. Its
endpoint will be the vector sum of those segments, ie at (x,y) = (sum(L_i
cos(A_i)), sum(L_i sin(A_i)).

To create a segmented path to a given (x,y), solve that expression for the
angles A_i. In R you can treat this as an optimisation problem - find a set
of angles A_i that minimise the distance of the end of the segmented path
from the target end point.

Here's some code that does that for a path from 0,0 to 0,1:

pathit <- function(segments){
    obj = function(angles){
        dxdy = dxdy(segments, angles)
        xerr = dxdy$dx-1
        yerr = dxdy$dy
        err = sqrt(xerr^2 + yerr^2)
        err
    }
    angles = runif(length(segments), 0 , 2*pi)
    optim(angles, obj)
}

dxdy = function(segments, angles){
    dx = sum(segments * cos(angles))
    dy = sum(segments * sin(angles))
    list(dx=dx, dy=dy)
}

plotsegs <- function(segments, angles){
    x = rep(NA, length(segments) +1)
    y = x
    x[1] = 0
    y[1] = 0
    for(i in 2:(length(segments)+1)){
        x[i] = x[i-1] + segments[i-1]*cos(angles[i-1])
        y[i] = y[i-1] + segments[i-1]*sin(angles[i-1])
    }
    cbind(x,y)
}

This is deliberately written naively for clarity.

To use, set the segment sizes, optimise, and then plot:

 s1 = c(.1,.3,.2,.1,.3,.3,.1)
 a1 = pathit(s1)
 plot(plotsegs(s1,a1$par),type="l")

which should show a path of seven segments from 0,0 to 0,1 - since the
initial starting values are random the model can find different solutions.
Run again for a different path.

To see what the space of paths looks like, do lots and overplot them:

  lots = lapply(1:1000, function(i)plotsegs(s1,pathit(s1)$par))
  plot(c(-.1,1.1),c(-1,1))
  p = lapply(lots, function(xy){lines(xy)})

this should show 1000 paths, and illustrates the "ellipse" of path
possibles that I mentioned in the previous email.

Sometimes the optimiser struggles to find a solution and so you should
probably test the output from optim for convergence and to make sure the
target function is close enough to zero for your purposes.  For the example
above most of the time the end point is about 1e-5  from (1,0) but for
harder problems such as s = rep(.1, 11) which only has 0.1 of extra "slack"
length, the error can be 0.02 and failed convergence. Possibly longer optim
runs would help or constraining the angles.

Anyway, interesting problem....

Barry







On Wed, Feb 6, 2019 at 8:23 PM Barry Rowlingson <b.rowlingson at gmail.com>
wrote:

  
  
6 days later
#
Hi everyone,

I wanted to generate random lines between two spatial points (random points
in polygons). The lines should consist of segments (9 segments), after
following suggestions I received to my previous post, I ended up using
trajectories. My aim now is to restrict these random trajectories (see
points in the attachment as an example) to the extent of a polygon (shown
in the attachment in yellow). So the trajectories should lie within the
polygon.
Please find the code I used below:
I first generated a data frame with my start and end point and used it to
generate a trajectory (straight line in attachment). I then generated a
random trajectory and shifted it so that the start and endpoints were equal
to my first trajectory. I transformed the matrices to a Spatial Data frame
to be able to plot it. Does anyone know how I can restrict the extent of
the trajectory that I randomly generate?

#coordinates of start and end point and time steps
x<-c(1649786,-1636500)
y<-c(-2902593,738500)
times<-c(0,9)
start_end_points<-as.data.frame(cbind(x,y,times))

# generate a trajectory between start and endpoint
trj<-TrajFromCoords(point_trj,xCol=1,yCol=2,timeCol=3, spatialUnits="m",
timeUnits="d")

# generate a random trajectory with 9 segments
track<-TrajGenerate(n=9, random=T, stepLength=TrajLength(trj)/9)

# To be able to use StartEndTranslate transform trajectories into matrices
ma1<-as.matrix(trj)
ma2<-as.matrix(track)

# Translate, rotate and scale the points of traj2 using traj1. The new traj
will have the same start and end points as traj1.
p<-StartEndTranslate(ma1,ma2)

# to check if this worked transform into SpatialDataFrame
# the start and end point of the new matrix are the same as the traj1 which
is good, the rest of the points changed as well
p_df<-as.data.frame(p)

# change lat / longs to nummeric
p_df$x<-as.numeric(p_df$x)
p_df$y<-as.numeric(p_df$y)

# generate a spatial dataframe to check if the StartEndTranslate function
worked by plotting it
xy<-p_df[,c(1,2)]

spdf<-SpatialPointsDataFrame(coords=xy, data=p_df, proj4string =
CRS(projection))

Thank you very much!

On Wed, Feb 6, 2019 at 5:42 PM Barry Rowlingson <b.rowlingson at gmail.com>
wrote:
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20190212/61f58f91/attachment.html>

-------------- next part --------------
A non-text attachment was scrubbed...
Name: restrict_trajectory_to_extent_of_polygon2.jpg
Type: image/jpeg
Size: 14874 bytes
Desc: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20190212/61f58f91/attachment.jpg>