Биннинг непрерывных значений с помощью round () создает артефакты

Question

В Python предположим, что у меня есть непрерывные переменные x и y, значения которых ограничены между 0 и 1 (чтобы было проще).Я всегда предполагал, что если я захочу преобразовать эти переменные в порядковые значения с ячейками, идущими как 0,0.01,0.02, ..., 0.98,0.99,1, можно просто округлить исходные значения до второй цифры.По какой-то причине, когда я это делаю, это оставляет артефакты.

Позвольте мне проиллюстрировать проблему (заметьте, однако, что мой вопрос не в том, как получить правильный график, а на самом деле, как сделать правильное объединение).Во-первых, это единственные модули, которые необходимы для воспроизведения проблемы:

import numpy as np
import matplotlib.pyplot as plt

Теперь предположим, что у нас есть непрерывные данные, сгенерированные следующим образом (другие процессы генерации данных также могут вызвать такую ​​же проблему):

# number of points drawn from Gaussian dists.:
n = 100000
x = np.random.normal(0, 2, n)
y = np.random.normal(4, 5, n)

# normalizing x and y to bound them between 0 and 1
# (it's way easier to illustrate the problem this way)
x = (x - min(x))/(max(x) - min(x))
y = (y - min(y))/(max(y) - min(y))

Затем давайте преобразуем x и y в порядковый номер в указанном выше интервале, просто применив некоторое округление.Затем давайте сохраним результаты в матрице x на y, чтобы построить ее тепловую карту для целей иллюстрации :

# matrix that will represent the bins. Notice that the
# desired bins are every 0.01, from 0 to 1, so 100 bins:
mtx = np.zeros([100,100])
for i in range(n):
    # my idea was that I could roughly get the bins by
    # simply rounding to the 2nd decimal point:
    posX = round(x[i], 2)
    posY = round(y[i], 2)
    mtx[int(posX*100)-1, int(posY*100)-1] += 1

Я ожидал, что вышеприведенное сработает, но когдаЯ строю содержимое матрицы mtx, на самом деле я получаю странные артефакты.Код:

# notice, however, the weird close-to-empty lines at
# 0.30 and 0.59 of both x and y. This happens regardless
# of how I generate x and y. Regardless of distributions
# or of number of points (even if it obviously becomes
# impossible to see if there are too few points):
plt.matshow(mtx, cmap=plt.cm.jet)
plt.show(block=False)

Дает мне:

Самое странное, что независимо от того, какой дистрибутив я использую для генерации x и y или какое семя я использую для ГСЧ, я всегда получаю одинаковые горизонтальные и вертикальные почти пустые линии в 0,30 и 0,59 как для x, так и y, довольно часто, причем линии, непосредственно параллельные этимпоказывает концентрацию точек (как вы видите на картинке).

Когда я печатаю значение по значению из этой матрицы в консоль, я могу фактически подтвердить, что те, которые соответствуют этим почти пустым строкам, действительно либо равны нулюили очень близко к нулю - в отличие от соседних точек.

Мой вопрос может быть более правильно разделен на 2 части:

1. Почему это произойдет?Я искренне хотел бы понять, что именно создает такую ​​проблему в этом простом коде.

2. Что было бы лучшим способом сгенерировать матрицу x с помощью y , который объединяет значения в соответствии с точками отсечения 0,0.01,0.02, ..., 0,98,0,99,1, не оставляя артефактов выше?

Если кто-то хочет легко захватить всеПример кода, использованного выше непосредственно в одной части, вот ссылка: https://www.codepile.net/pile/VLAq4kLp

ПРИМЕЧАНИЕ: я не хочу найти правильный способ построения.Я хочу найти myeself правильный способ генерации "матрицы значений", представленной на рисунке выше.Я знаю, что есть и другие способы составления графика тепловой карты без артефактов, например, с использованием plt.matshow(mtx, cmap=plt.cm.jet); plt.show(block=False) или plt.hist2d(x, y, bins=100).Что я спрашиваю, так это где проблема в самой генерации моей матрицы, которая создает эти почти нулевые элементы.

Answer 1

Проблема с вашим методом - ошибка с плавающей запятой.Это становится очевидным, когда вы пытаетесь превратить ваше округленное число в целое число.Рассмотрим следующую функцию (по сути, то, что вы делаете с каждым из ваших случайных чисел):

def int_round(a):
     r = round(a, 2)
     rh = r*100
     i = int(rh)
     print(r, rh, i)


int_round(0.27)
#prints: 0.27 27.0 27

int_round(0.28)
#prints: 0.28 28.000000000000004 28

int_round(0.29)
#prints: 0.29 28.999999999999996 28

int_round(0.30)
#prints: 0.3 30.0 30

Как видите, из-за ошибки с плавающей запятой после округления 0,28 и 0,29 и умножения на 100,и 0.28, и 0.29 заканчиваются целым числом 28.(Это потому, что int() всегда округляется, поэтому 28.99999999999 становится 28).

Решением может быть округление значения после умножения на 100:

def round_int(a):
    ah = a*100
    rh = round(ah, 2)
    i = int(rh)
    print(ah, rh, i)

round_int(0.27)
#prints: 27.0 27.0 27

round_int(0.28)
#prints: 28.000000000000004 28.0 28

round_int(0.29)
#prints: 28.999999999999996 29.0 29

round_int(0.30)
#prints: 30.0 30.0 30

Обратите внимание, что в этом случае0.29 исправлено и преобразовано в 29.

Применение этой логики к вашему коду: мы можем изменить цикл for на:

mtx = np.zeros([101, 101])

for i in range(n):
    # my idea was that I could roughly get the bins by
    # simply rounding to the 2nd decimal point:
    posX = np.round(100*x[i], 2)
    posY = np.round(100*y[i], 2)
    mtx[int(posX), int(posY)] += 1

Обратите внимание на увеличение количества бинов до101 для учета последнего бина, когда x = 1 или y = 1.Кроме того, здесь вы можете видеть, что при умножении x[i] и y[i] на 100 перед округлением биннинг происходит правильно:

Answer 2

Я не знаю, как точно ответить на ваш первый вопрос.Но для хранения предметов я также использую pandas.cut .Для вашего решения вы можете сделать

import pandas as pd
bins = [v / 100. for v in range(100)
bucketed = pd.cut(x, bins)

bucketed, а затем указать, к какому интервалу относится каждая точка данных

. Для справки вот достойное руководство по этому вопросу http://benalexkeen.com/bucketing-continuous-variables-in-pandas/

Answer 3

Проблема может быть легко решена с помощью np.histogram2d(x,y, bins=100).

. Остаток ответа - показать, где ручные алгоритмы не работают:

Считайте, что численно

0.56*100 == 56.00000000000001    -> int(0.56*100) == 56
0.57*100 == 56.99999999999999    -> int(0.57*100) == 56
0.58*100 == 57.99999999999999    -> int(0.58*100) == 57
0.59*100 == 59.00000000000000    -> int(0.59*100) == 59

так, что число 58 просто не будет присутствовать в вашей индексации, а число 56 будет появляться в два раза чаще (для равномерного распределения).

Вместо этого вы можете сначала умножить, а затем усечь до целого числа.Также обратите внимание, что последний лоток должен быть закрыт, так что значение 1 добавляется в лоток с индексом 99.

mtx = np.zeros([100,100])
for i in range(n):
    posX = int(x[i]*100)
    posY = int(y[i]*100)
    if posX == 100:
        posX = 99
    if posY == 100:
        posY = 99
    mtx[posX, posY] += 1

Это будет определять лотки через ребра, т. Е. Первый интервал варьируется отОт 0 до 1 и т. Д. При вызове imshow / matshow вам нужно будет принять это во внимание, установив экстент.

plt.matshow(mtx, cmap=plt.cm.jet, extent=(0,100,0,100))

Answer 4

На данный момент я могу только правильно ответить на ваш второй вопрос, поскольку я все еще ищу ошибку в первой части.

Итак, вот стандартное решение, которое вы бы выбрали для такого типа, как вы.хочу (предполагая, что x и y, которые вы упомянули ранее):

h = plt.hist2d(x, y, bins=100)

, дающий

, который являетсяСетка 100x100.

Переменная h теперь содержит нужную матрицу, а также ячейки, найденные matplotlib.plt.matshow(h[0]) показывает ту же матрицу, что и на рисунке, которую возвращает matplotlib.Как уже упоминалось в комментариях: вы можете получить те же результаты (но без автоматического графика), вызвав

h = np.histogram2d(x, y, bins=100)

Тем не менее, ваш алгоритм не может быть верным, потому что вы фактически подсчитываете количество элементов по краям, а не между их, поэтому вы получаете 101 предмет в каждом направлении.Вы можете увидеть проблему, например, posX==0: Тогда int(posX*100)-1 дает -1.