Calculate a measure of "almost multiplicity" from a list of noisy values

Question

Question

Suppose I have a list of values ​​that have a total number greater than 1. For example, take a numerical number equal to 3 and form a set of multiple values ​​of this value:

harmonicList = [3,6,6,3,3,9,27,3,15,18,9]

Now I add some noise:

harmonicList = [ v + (random() * 2 - 1 ) * 0.1 for v in harmonicList]

I am looking for an algorithm that will return a number around 1.0 when the elements in the list are close to a multiple of the total value, but close to 0.0 when the numbers are not close to a multiple - for example, for example, when the list is a set of primes.

Is there such a measure of "almost multiplicity"?

Why I want to solve this problem

I'm currently trying to detect a Chessboard in a screenshot using Hough Transform. Sometimes the case is perfect and it works very well:

But sometimes not:

, . , , ( , ). , , "", .

, , , . , , , , .

Answer 1

Python . , score - , - (x ). 0 1, :

import numpy as np
import scipy.stats as stats
import matplotlib.pyplot as plt


primes = np.array([2,  3,  5,  7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43,
                   47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97])

def harmonic(divisor, size=10):
    return np.random.randint(1, 10, size=size) * divisor

def prime_sample(size=10):
    return np.random.choice(primes, size=size)

def noisy(x, amount=0.1):
    return x + (np.random.random(size=len(x)) * 2 - 1) * amount

def prob(x, mean, sd):
    return stats.norm.pdf(x, loc=mean, scale=sd)

def score(x, multiplier, offset, kmax=20):
    k = np.arange(kmax)
    means = (k * multiplier + offset)[:, None]
    closest_multiple = (np.abs(x - means).argmin(axis=0)) * multiplier
    result = np.exp(-((x - closest_multiple)**2).sum())
    return result

def fit(x, multipliers, offsets, kmax=20, sd=0.2):
    "sd is the standard deviation of the noise"
    k = np.arange(kmax)
    M, O, K = np.meshgrid(multipliers, offsets, k, indexing='ij')
    means = (K * M + O)[..., None]
    p = prob(x, means, sd)
    # sum over the K axis, take the log, sum over x axis
    L = np.log(p.sum(axis=-2)).sum(axis=-1)
    # find the location of maximum log likelihood
    i, j = np.unravel_index(L.argmax(), L.shape)
    max_L = L[i, j]
    multiplier = multipliers[i]
    offset = offsets[j]
    return dict(loglikelihood=L, max_L=max_L,
                multiplier=multiplier, offset=offset,
                score=score(x, multiplier, offset, kmax))

multipliers = np.linspace(3, 10, 100)
offsets = np.linspace(-1.5, 1.5, 50)
X, Y = np.meshgrid(multipliers, offsets, indexing='ij')
tests = [([12,  8, 28, 20, 32, 12, 28, 16,  4, 12], 1),
         ([3, 5, 7, 11, 13, 27, 54, 57], 0),
         (noisy(harmonic(3, size=20)), 1),
         (noisy(prime_sample()), 0)]

for x, expected in tests:
    result = fit(x, multipliers, offsets, kmax=20)
    Z = result['loglikelihood']
    plt.contourf(X, Y, Z)
    plt.xlabel('multiplier')
    plt.ylabel('offset')
    plt.scatter(result['multiplier'], result['offset'], s=20, c='red')
    plt.title('score = {:g}, expected = {:g}'
              .format(result['score'], expected))
    plt.show()

x = [12, 8, 28, 20, 32, 12, 28, 16, 4, 12]:

x = [3, 5, 7, 11, 13, 27, 54, 57]:

Answer 2

, - , :

for i in range(1,int(round(min(harmonicList)))+1): #makes no sense to look at bigger numbers
    temp_sum = 0
    for item in harmonicList:
        item = item % i
        temp_sum += min(item, i-item) #because noise can be either plus or minus
    if i == 1:
        comparator = temp_sum
    elif temp_sum == comparator:
        print(i, "is a common denominator")
#prints "3 is a common denominator"

:

primeList = [1,3,5,7,11,13,17]

for , .

for _ in range(1000):
    primeList = [1,3,5,7,11,13,17]
    primeList = [ v + (random() * 2 - 1 ) * 0.1 for v in primeList]
    for i in range(1,int(round(min(primeList)))+1):
        temp_sum = 0
        for item in primeList:
            item = item % i
            temp_sum += min(item, i-item)
        if i == 1:
            comparator = temp_sum
        elif temp_sum == comparator:
            print(i, "is a common denominator")
#No prints

Answer 3

Not very clean (hopefully someone improves), but this does the job:

from operator import itemgetter
def near_multiplicity(num_list):
    rounded_num_list = [round(num) for num in num_list]
    factors = {}
    for num in rounded_num_list:
        for factor in range(2, int(num)+1):
            if int(num) % factor == 0:
                factors[factor] = factors.get(factor, 0) + 1
    sorted_factors = sorted(factors.items(), key=itemgetter(1), reverse=True)
    if sorted_factors[0][1] == 1:
        return 0
    best_factor = sorted_factors[0][0]
    noise = 0
    for num in num_list:
        distortion = num % best_factor
        noise += min(distortion, best_factor - distortion)
    average_noise = noise / len(num_list)
    return 1 - average_noise