Search code examples
pythoncluster-analysisspatialdistribution

Is there an algorithm for generating random samples in a dispersed distribution, as if they were actively repelling one another?


I am attempting to create a script for generating randomly distributed samples with a dispersed distribution - as if they are actively repelling one another - such as with the example in the linked image.

Example of spatial point patterns

That image pertains to Ripley's K algorithm, which is useful for analysing the nature of an existing sample distribution - but I'm attempting to find an algorithm to actually generate points in a dispersed distribution, and I'm not having any luck.

I attempted to write my own algorithm for this purpose, but it's very ineffective and crudely put together. This algorithm just randomly picks points in space and checks whether said point is within N units of any existing sample, with N gradually reducing over time to avoid circumstances where the code would hang. At high amounts of samples, the distributions generated from this code are no better than random.

    for n in range(numSamples):
        tempIdealDistance = idealDistance #Units of distance we're attempting to separate each sample by
        newPosition = []
        searchingForLocation = True
        while searchingForLocation:
            #Looking for a location that doesn't overlap with any existing samples
            newPosition = [random.randint(lowBounds,highBounds),random.randint(lowBounds,highBounds),random.randint(lowBounds,highBounds)]
            withinRange = False
            for s in samples:
                #Iterate through the existing list of samples, looking for overlaps
                existingSamplePosition = np.array(s.GetPos())
                newSamplePosition = np.array(newPosition)
                if np.linalg.norm(existingSamplePosition - newSamplePosition) < tempIdealDistance:
                    #If we overlap with any existing samples, this is not a valid position
                    withinRange = True
                    break
            if withinRange == False:
                #If no samples are within range of our new one, it fits the regular distribution
                searchingForLocation = False
                samp = Point(newPosition[0],newPosition[1],newPosition[2])
                samp.SetRatio(GetRatio(TiSrB,TiSrA,samp.GetPos(),position),randomNoise)
                samples.append(samp)
            else:
                #To avoid infinite loops, we gradually decrease the search distance
                tempIdealDistance = tempIdealDistance - 1

Could anyone help me find a more appropriate algorithm to what I'm looking for?


Solution

  • I could think about couple of options:

    1. Poisson disc sampling. See https://www.jasondavies.com/poisson-disc/ for details. Basically samples would be separated by minimum distance.

    2. Quasi-random (e.g. Sobol) sequences. There is a theorem for Sobol sequences which states that minimum distance would be sqrt(d)/2 1/N, where d is problem dimensions, and N is number of samples. See https://en.wikipedia.org/wiki/Sobol_sequence for details