Re: [R] quasi-random sequences
Hi again,
I've had a go at Prof Ripley's suggestion (Strauss process, code
below). It works great for my limited purpose (qualitative drawing,
really), I can just add a few mild concerns,
- ideally the hard core would not be a fixed number, but a
distribution of sizes (ellipses).
- I could not quite understand the link between the window size,
intensity, and number of elements returned. Is there a simple
relation I've missed?
- in general, I just have no clue how rStrauss works, I guess Prof
Ripley's first reference cited in ?Strauss would be useful for that
matter as I could not find anything with google but a C code looking
a bit alike some genetic algorithm (death and birth of objects).
Am I on the right track in thinking that the idea is not so
dissimilar with the following approach:
1) create a first random distribution (Poisson, apparently)
2) give birth and death to new objects
3) carry on with this evolution until the sample distribution
respects the given intensity and interaction parameters (e.g, for
hard core, no overlap is permitted within the radius)
This looks a bit like the idea of an electrostatic problem,
1) start with an initial configuration
2) define the total potential energy as a cost function to be
minimized, sum of a boundary (repulsive term), and short distance
repulsion
3) optimize the positions with respect to this objective function
I tried unsuccessfully this second approach (see below), but I'm sort
of concerned about the handling of all these parameters by optim().
rStrauss works beautifully anyway !
Thanks again,
baptiste
##
# Code
##
library(spatstat) # Strauss process
## ellipse function obtained from R mailing list ##
## could use library(car) as an alternative
ellipse - function(x,y,width,height,theta, npoints=100,plot=F)
{
# x = x coordinate of center
# y = y coordinate of center
# width = length of major axis
# height = length of minor axis
# theta = rotation
# npoints = number of points to send to polygon
# plot = if TRUE, add to current device
# = if FALSE, return list of components
a - width/2
b - height/2
theta-theta*pi/180
xcoord - seq(-a,a,length=npoints)
ycoord.neg - sqrt(b^2*(1-(xcoord)^2/a^2))
ycoord.pos - -sqrt(b^2*(1-(xcoord)^2/a^2))
xx - c(xcoord,xcoord[npoints:1])
yy - c(ycoord.neg,ycoord.pos)
x.theta - xx*cos(2*pi-theta)+yy*sin(2*pi-theta)+x
y.theta - yy*cos(2*pi-theta)-xx*sin(2*pi-theta)+y
if(plot)
invisible(polygon(x.theta,y.theta,col=black))
else
invisible(list(coords=data.frame(x=x.theta,y=y.theta),
center=c(x,y),
theta=theta))
}
getEllipse-function(dataf,plot=TRUE){
ellipse(dataf[1],dataf[2],dataf[3],dataf[4],dataf[5],plot=plot)
}
N - 200 # initial number (loosely related to the actual number of
elements)
##
## positions from Strauss hard core process ##
##
positions - rStrauss(beta = 0.1, gamma = 0, R = 2, W= square(N/2))
# beta: intensity
# gamma: interaction. 0 : hard core (no overlap), 1: complete
randomness
# R: radius of interaction (size of the core)
# W: window
N2-length(positions$x)# number of elements
rand.angles - rnorm(N2,sd=45,mean=45) # ellipse angles
rand.a - rnorm(N2,sd=0.4,mean=1) # ellipse semi-axes
rand.b - rnorm(N2,sd=0.4,mean=1) #
par(bty=n)
plot(cbind(positions$x,positions$y) ,cex = rand.a,axes=F,
xlab=,ylab=,t=n) # just circles: type=p
ell.pos-cbind(positions$x,positions$y,rand.a,rand.b,rand.angles)
apply(ell.pos,1,getEllipse) - b.quiet # draw ellipses
##
##
##
## optimizing an electrostatic potential problem (not working)
##
N -200
p0 - matrix(rnorm(2*N),ncol=2) # random starting point
# boundary repulsion potential #
x - seq(-1.2,1.2,l=N)
delta - 0.1
# 1D example #
y - -1*( 1 - exp(-abs(x)/delta) - exp(abs(x)/delta))
#plot(x,y,t=l)
# 2D example
xy - expand.grid(x = x, y=x)
z.x - -1*( 1 - exp(-abs(xy[,1])/delta) - exp(abs(xy[,1])/delta))
z.y - -1*( 1 - exp(-abs(xy[,2])/delta) - exp(abs(xy[,2])/delta))
z - matrix(z.x+z.y, ncol=length(x))
#image(x=x,y=x,z=z)
boundary - function(xy = p0[1,]){
z.x - -1*( 1 - exp(-abs(xy[1])/delta) - exp(abs(xy[1])/delta))
z.y - -1*( 1 - exp(-abs(xy[2])/delta) - exp(abs(xy[2])/delta))
z.x + z.y
}
dist.2d - function(x = c(1, 0 ),y = c(0,1) , w = 0.1){
(drop((x - y) %*% (x - y)) + w^2 )^(-3/2)
}
obj - function(p = p0){
bound - sum(apply(p,1,boundary)) # boundary term
rel.positions - expand.grid(x=p,y=p) # all relative positions
bound + sum(apply(rel.positions,1,dist.2d)) #
Re: [R] quasi-random sequences
baptiste Auguié wrote: Dear list useRs, I have to generate a random set of coordinates (x,y) in [-1 ; 1]^2 for say, N points. At each of these points is drawn a circle (later on, an ellipse) of random size, as in: The quasi-random sequences are useful for integration, but they're not random in the sense of being unpredictable. Here's an idea for a sequence that might work for you: divide the square into an n x n grid of smaller squares. Randomly select a smaller square, then uniformly select a point within it. (Or just permute the list of all n^2 squares, and proceed through the permuted list.) This won't guarantee a separation, but it will tend to lead to one. If you want a guarantee, then only select from squares with even coordinates. Once you've selected all the squares in the grid, go to a 2n x 2n grid, and repeat. (If you are doing the even coordinates modification, you'll have had to select non-uniformly for this to work.) Duncan Murdoch N - 100 positions - matrix(rnorm(2 * N, mean = 0 , sd= 0.5), nrow=N) sizes-rnorm(N, mean = 0 , sd= 1) plot(positions,type=p,cex=sizes) My problem is to avoid collisions (overlap, really) between the points. I would like some random pattern, but with a minimum exclusion distance. In looking up Numerical recipes in C, I found out about some Sobol quasi-random sequences, which one can call from the gsl package, library(gsl) g - qrng_alloc(type=sobol,dim=2) qrng_get(g,n= N) -xy plot((xy),t=p,cex=0.5) but this does not look very random: I clearly see some pattern (diagonals, etc...), and even the non-overlapping condition is not impressive. One (painful) way I can foresee is to check the distance between each symbol and the others, and move the overlapping ones in a recursive manner. Before delving into this, I wanted to check I'm not overlooking something in the rgl quasi-random sequences, or missing a more obvious way to generate such patterns. Perhaps solving an electrostatic problem with a potential both attractive at long distances and repulsive at short distances is a better way? I have a vague recollection of hearing that somewhere to position points evenly on a sphere. Thanks for any comment / suggestion, Baptiste _ Baptiste Auguié Physics Department University of Exeter Stocker Road, Exeter, Devon, EX4 4QL, UK Phone: +44 1392 264187 http://newton.ex.ac.uk/research/emag http://projects.ex.ac.uk/atto __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] quasi-random sequences
On Sat, 26 Apr 2008, baptiste Auguié wrote: Dear list useRs, I have to generate a random set of coordinates (x,y) in [-1 ; 1]^2 for say, N points. At each of these points is drawn a circle (later on, an ellipse) of random size, as in: N - 100 positions - matrix(rnorm(2 * N, mean = 0 , sd= 0.5), nrow=N) sizes-rnorm(N, mean = 0 , sd= 1) plot(positions,type=p,cex=sizes) My problem is to avoid collisions (overlap, really) between the points. I would like some random pattern, but with a minimum exclusion distance. In looking up Numerical recipes in C, I found out about some Sobol quasi-random sequences, which one can call from the gsl package, That is a hard-core point process, e.g. a special case of a Strauss process. You will find ways to simulate such a process (there are several processes, and several ways for most) in the various spatial packages, including 'spatial' itself. I think you have misunderstood the point of quasi-random sequences (which, given the exposition in the edition of NR I have, would be easy to do). library(gsl) g - qrng_alloc(type=sobol,dim=2) qrng_get(g,n= N) -xy plot((xy),t=p,cex=0.5) but this does not look very random: I clearly see some pattern (diagonals, etc...), and even the non-overlapping condition is not impressive. One (painful) way I can foresee is to check the distance between each symbol and the others, and move the overlapping ones in a recursive manner. Before delving into this, I wanted to check I'm not overlooking something in the rgl quasi-random sequences, or missing a more obvious way to generate such patterns. Perhaps solving an electrostatic problem with a potential both attractive at long distances and repulsive at short distances is a better way? I have a vague recollection of hearing that somewhere to position points evenly on a sphere. Thanks for any comment / suggestion, Baptiste _ Baptiste Auguié Physics Department University of Exeter Stocker Road, Exeter, Devon, EX4 4QL, UK Phone: +44 1392 264187 http://newton.ex.ac.uk/research/emag http://projects.ex.ac.uk/atto __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595__ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] quasi-random sequences
You seem to have ambiguous requirements. First, you want equidistribution for a packing structure, which would suggest closest packing or quasirandom sequences, as you have tried. But then you are disturbed by the packing structure, because it gives a perceivable pattern, so you wish to randomize it, but under some unspecified condition of equidistribution (your electrostatic repulsion algorithm). You obviously have some constraints on selection you have not quantified yet. E.g., circles of unspecified radii cannot overlap. You should also realize that any distribution of centers under such constraints will always exhibit structure due to the constraints. Is your problem simply to give an appearance of randomness to the casual observer, or something more definite? You also need to say something about the packing density involved. On the face of it, you with At 06:22 AM 4/26/2008, baptiste Auguié wrote: Dear list useRs, I have to generate a random set of coordinates (x,y) in [-1 ; 1]^2 for say, N points. At each of these points is drawn a circle (later on, an ellipse) of random size, as in: N - 100 positions - matrix(rnorm(2 * N, mean = 0 , sd= 0.5), nrow=N) sizes-rnorm(N, mean = 0 , sd= 1) plot(positions,type=p,cex=sizes) My problem is to avoid collisions (overlap, really) between the points. I would like some random pattern, but with a minimum exclusion distance. In looking up Numerical recipes in C, I found out about some Sobol quasi-random sequences, which one can call from the gsl package, library(gsl) g - qrng_alloc(type=sobol,dim=2) qrng_get(g,n= N) -xy plot((xy),t=p,cex=0.5) but this does not look very random: I clearly see some pattern (diagonals, etc...), and even the non-overlapping condition is not impressive. One (painful) way I can foresee is to check the distance between each symbol and the others, and move the overlapping ones in a recursive manner. Before delving into this, I wanted to check I'm not overlooking something in the rgl quasi-random sequences, or missing a more obvious way to generate such patterns. Perhaps solving an electrostatic problem with a potential both attractive at long distances and repulsive at short distances is a better way? I have a vague recollection of hearing that somewhere to position points evenly on a sphere. Thanks for any comment / suggestion, Baptiste _ Baptiste Auguié Physics Department University of Exeter Stocker Road, Exeter, Devon, EX4 4QL, UK Phone: +44 1392 264187 http://newton.ex.ac.uk/research/emag http://projects.ex.ac.uk/atto __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Robert A. LaBudde, PhD, PAS, Dpl. ACAFS e-mail: [EMAIL PROTECTED] Least Cost Formulations, Ltd.URL: http://lcfltd.com/ 824 Timberlake Drive Tel: 757-467-0954 Virginia Beach, VA 23464-3239Fax: 757-467-2947 Vere scire est per causas scire __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] quasi-random sequences
baptiste Auguié ba208 at exeter.ac.uk writes: Dear list useRs, You might be interested to apply the Hammersley or Halton point sets that are often used in numerical integration or Differential Evolution. These pseudo-random distributions are both uniform and irregular, but have a kind of minimum resolution There is an implementation of Halton Sequences in the often overlooked 'sfsmisc' package, see the 'QUnif()' function there. The help includes a visualization example dispaying the behavior I think you were looking for. Hans Werner I have to generate a random set of coordinates (x,y) in [-1 ; 1]^2 for say, N points. At each of these points is drawn a circle (later on, an ellipse) of random size, [...] My problem is to avoid collisions (overlap, really) between the points. I would like some random pattern, but with a minimum exclusion distance. In looking up Numerical recipes in C, I found out about some Sobol quasi-random sequences, which one can call from the gsl package, [...] but this does not look very random: I clearly see some pattern (diagonals, etc...), and even the non-overlapping condition is not impressive. One (painful) way I can foresee is to check the distance between each symbol and the others, and move the overlapping ones in a recursive manner. Before delving into this, I wanted to check I'm not overlooking something in the rgl quasi-random sequences, or missing a more obvious way to generate such patterns. Perhaps solving an electrostatic problem with a potential both attractive at long distances and repulsive at short distances is a better way? I have a vague recollection of hearing that somewhere to position points evenly on a sphere. Thanks for any comment / suggestion, Baptiste __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] quasi-random sequences
Thank you all for the great suggestions and comments. As two of you pointed out, the problem was not well defined (who said a well-posed problem is a problem solved?), and also it seems to be a very wide topic. I've had an interesting reading discussing the similarities between half-toning in black and white printing, and the quasi Monte Carlo integration technique [1]. I'm trying the suggested ways plus a few other things, I'll put together a summary of my discoveries when I've clarified it all. Thanks again, baptiste [1] Halftoning and Quasi-Monte Carlo, K. M. Hanson, Sensitivity Analysis of Model Output, K. M. Hanson and F. M. Hemez, eds., pp. 430-442 (Los Alamos Research Library, 2005) _ Baptiste Auguié Physics Department University of Exeter Stocker Road, Exeter, Devon, EX4 4QL, UK Phone: +44 1392 264187 http://newton.ex.ac.uk/research/emag http://projects.ex.ac.uk/atto __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] quasi-random sequences
You might want to shuffle coordinates independently to get rid of the diagonals. Otherwise what quasi-random sequence guarantee are upper boundaries on the coverage errors, but not anything nice-looking and irregular. Sobol' sequences, even though they are theoretically superior to some others (e.g., Halton sequences more popular among economists), are especially nasty in producing bands and bricks on the low dimensional plots. Among statisticians, Art Owen from Stanford is almost the only one interested in this sort of stuff (referred to as quasi-Monte Carlo, in his field(s)). You might have better luck on a physics list with a question like yours. On Sat, Apr 26, 2008 at 5:22 AM, baptiste Auguié [EMAIL PROTECTED] wrote: Dear list useRs, I have to generate a random set of coordinates (x,y) in [-1 ; 1]^2 for say, N points. At each of these points is drawn a circle (later on, an ellipse) of random size, as in: N - 100 positions - matrix(rnorm(2 * N, mean = 0 , sd= 0.5), nrow=N) sizes-rnorm(N, mean = 0 , sd= 1) plot(positions,type=p,cex=sizes) My problem is to avoid collisions (overlap, really) between the points. I would like some random pattern, but with a minimum exclusion distance. In looking up Numerical recipes in C, I found out about some Sobol quasi-random sequences, which one can call from the gsl package, library(gsl) g - qrng_alloc(type=sobol,dim=2) qrng_get(g,n= N) -xy plot((xy),t=p,cex=0.5) but this does not look very random: I clearly see some pattern (diagonals, etc...), and even the non-overlapping condition is not impressive. One (painful) way I can foresee is to check the distance between each symbol and the others, and move the overlapping ones in a recursive manner. Before delving into this, I wanted to check I'm not overlooking something in the rgl quasi-random sequences, or missing a more obvious way to generate such patterns. Perhaps solving an electrostatic problem with a potential both attractive at long distances and repulsive at short distances is a better way? I have a vague recollection of hearing that somewhere to position points evenly on a sphere. -- Stas Kolenikov, also found at http://stas.kolenikov.name Small print: I don't check Gmail account regularly. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
