python中的圆形/极坐标直方图

Posted 2023-02-14

【中文标题】python中的圆形/极坐标直方图

【英文标题】：Circular / polar histogram in python

【发布时间】：2014-04-29 00:47:44

【问题描述】：

我有周期性数据，它的分布最好围绕一个圆圈可视化。现在的问题是如何使用matplotlib 进行可视化？如果没有，用 Python 可以轻松完成吗？

这里我生成了一些我想用圆形直方图可视化的示例数据：

import matplotlib.pyplot as plt
import numpy as np

# Generating random data
a = np.random.uniform(low=0, high=2*np.pi, size=50)

Mathematica 的 SX 问题中有几个示例。

我想生成一个类似于以下内容之一的图：

【问题讨论】：

我没有关注...我必须证明我是从头开始编写的，还是应该要求人们从头开始编写？

这会让你开始：matplotlib.org/examples/pie_and_polar_charts/…

@PaulH，非常感谢 :)

【参考方案1】：

以图库中的this 示例为基础，您可以这样做

import numpy as np
import matplotlib.pyplot as plt

N = 80
bottom = 8
max_height = 4

theta = np.linspace(0.0, 2 * np.pi, N, endpoint=False)
radii = max_height*np.random.rand(N)
width = (2*np.pi) / N

ax = plt.subplot(111, polar=True)
bars = ax.bar(theta, radii, width=width, bottom=bottom)

# Use custom colors and opacity
for r, bar in zip(radii, bars):
    bar.set_facecolor(plt.cm.jet(r / 10.))
    bar.set_alpha(0.8)

plt.show()

当然，有很多变体和 tweeks，但这应该可以帮助您入门。

一般来说，浏览matplotlib gallery 通常是一个不错的起点。

在这里，我使用bottom 关键字将中心留空，因为我想我看到了您之前的问题，其图表更像我所拥有的，所以我认为这就是您想要的。要获得上面显示的完整楔形，只需使用 bottom=0（或将其省略，因为 0 是默认值）。

【讨论】：

你知道怎么开始左边的0度而不是180度吗？

我认为ax.set_theta_zero_location("W")。 （不过，一般来说，最好提出一个新问题而不是发表评论。这样，后续、更改、示例图等都可以添加。）

非常感谢，它成功了，虽然它使底部的 90 度和顶部的 180 度。

啊我用ax.set_theta_direction(-1)！

ax.set_theta_offset(offset_in_radians) 更改 matplotlib 2.1.0 中的方向

【参考方案2】：

快速解答

使用下面我写的函数circular_hist()。

默认情况下，此函数绘制与面积成比例的频率，而不是半径（此决定背后的原因在下面的“更长的形式答案”下提供）。

def circular_hist(ax, x, bins=16, density=True, offset=0, gaps=True):
    """
    Produce a circular histogram of angles on ax.

    Parameters
    ----------
    ax : matplotlib.axes._subplots.PolarAxesSubplot
        axis instance created with subplot_kw=dict(projection='polar').

    x : array
        Angles to plot, expected in units of radians.

    bins : int, optional
        Defines the number of equal-width bins in the range. The default is 16.

    density : bool, optional
        If True plot frequency proportional to area. If False plot frequency
        proportional to radius. The default is True.

    offset : float, optional
        Sets the offset for the location of the 0 direction in units of
        radians. The default is 0.

    gaps : bool, optional
        Whether to allow gaps between bins. When gaps = False the bins are
        forced to partition the entire [-pi, pi] range. The default is True.

    Returns
    -------
    n : array or list of arrays
        The number of values in each bin.

    bins : array
        The edges of the bins.

    patches : `.BarContainer` or list of a single `.Polygon`
        Container of individual artists used to create the histogram
        or list of such containers if there are multiple input datasets.
    """
    # Wrap angles to [-pi, pi)
    x = (x+np.pi) % (2*np.pi) - np.pi

    # Force bins to partition entire circle
    if not gaps:
        bins = np.linspace(-np.pi, np.pi, num=bins+1)

    # Bin data and record counts
    n, bins = np.histogram(x, bins=bins)

    # Compute width of each bin
    widths = np.diff(bins)

    # By default plot frequency proportional to area
    if density:
        # Area to assign each bin
        area = n / x.size
        # Calculate corresponding bin radius
        radius = (area/np.pi) ** .5
    # Otherwise plot frequency proportional to radius
    else:
        radius = n

    # Plot data on ax
    patches = ax.bar(bins[:-1], radius, zorder=1, align='edge', width=widths,
                     edgecolor='C0', fill=False, linewidth=1)

    # Set the direction of the zero angle
    ax.set_theta_offset(offset)

    # Remove ylabels for area plots (they are mostly obstructive)
    if density:
        ax.set_yticks([])

    return n, bins, patches

示例用法：

import matplotlib.pyplot as plt
import numpy as np

angles0 = np.random.normal(loc=0, scale=1, size=10000)
angles1 = np.random.uniform(0, 2*np.pi, size=1000)

# Construct figure and axis to plot on
fig, ax = plt.subplots(1, 2, subplot_kw=dict(projection='polar'))

# Visualise by area of bins
circular_hist(ax[0], angles0)
# Visualise by radius of bins
circular_hist(ax[1], angles1, offset=np.pi/2, density=False)

更长的答案

我总是建议在使用圆形直方图时要小心，因为它们很容易误导读者。

特别是，我建议远离按比例绘制 频率 和 半径 的圆形直方图。我推荐这个，因为思维会受到垃圾箱区域的极大影响，而不仅仅是它们的径向范围。这类似于我们习惯于解释饼图的方式：按区域。

因此，我建议不要使用 bin 的 radial 范围来可视化它包含的数据点数量，而是按区域可视化点数。

问题

考虑将给定直方图 bin 中的数据点数量加倍的后果。在频率和半径成比例的圆形直方图中，这个 bin 的半径将增加 2 倍（因为点数增加了一倍）。但是，这个垃圾箱的面积将增加 4 倍！这是因为 bin 的面积与半径的平方成正比。

如果这听起来还不是什么大问题，让我们用图形来看看：

以上两个图都可视化了相同的数据点。

在左侧图中，很容易看出 (0, pi/4) bin 中的数据点数量是 (-pi/4, 0) bin 中的两倍。

但是，请看一下右侧图（频率与半径成正比）。乍一看，您的思想受到垃圾箱面积的极大影响。如果您认为 (0, pi/4) bin 中的点数是 (-pi/4, 0) bin 中的点数的 多于 倍，这是可以原谅的。然而，你会被误导。只有在仔细检查图形（和径向轴）时，您才会意识到 (0, pi/4) bin 中的数据点是 (-pi) 的两倍/4, 0) 箱。 不超过图表最初建议的两倍。

上面的图形可以用下面的代码重新创建：

import numpy as np
import matplotlib.pyplot as plt
plt.style.use('seaborn')

# Generate data with twice as many points in (0, np.pi/4) than (-np.pi/4, 0)
angles = np.hstack([np.random.uniform(0, np.pi/4, size=100),
                    np.random.uniform(-np.pi/4, 0, size=50)])

bins = 2

fig = plt.figure()
ax = fig.add_subplot(1, 2, 1)
polar_ax = fig.add_subplot(1, 2, 2, projection="polar")

# Plot "standard" histogram
ax.hist(angles, bins=bins)
# Fiddle with labels and limits
ax.set_xlim([-np.pi/4, np.pi/4])
ax.set_xticks([-np.pi/4, 0, np.pi/4])
ax.set_xticklabels([r'$-\pi/4$', r'$0$', r'$\pi/4$'])

# bin data for our polar histogram
count, bin = np.histogram(angles, bins=bins)
# Plot polar histogram
polar_ax.bar(bin[:-1], count, align='edge', color='C0')

# Fiddle with labels and limits
polar_ax.set_xticks([0, np.pi/4, 2*np.pi - np.pi/4])
polar_ax.set_xticklabels([r'$0$', r'$\pi/4$', r'$-\pi/4$'])
polar_ax.set_rlabel_position(90)

解决方案

由于我们受到圆形直方图中 bin 的 面积 的极大影响，我发现确保每个 bin 的面积与其中的观察次数成正比会更有效，而不是的半径。这类似于我们习惯于解释饼图的方式，其中面积是感兴趣的数量。

让我们使用我们在上一个示例中使用的数据集来根据面积而不是半径来重现图形：

我相信读者在第一眼看到这张图片时被误导的可能性较小。

但是，在绘制面积与半径成正比的圆形直方图时，我们的缺点是您永远不会知道 (0, pi/4) 中有 恰好 两倍的点bin 比 (-pi/4, 0) bin 仅仅通过观察这些区域。虽然，您可以通过用相应的密度注释每个 bin 来解决这个问题。我认为这个缺点比误导读者更可取。

当然，我会确保在此图旁边放置一个信息丰富的标题，以解释此处我们将频率可视化为面积，而不是半径。

以上地块创建为：

fig = plt.figure()
ax = fig.add_subplot(1, 2, 1)
polar_ax = fig.add_subplot(1, 2, 2, projection="polar")

# Plot "standard" histogram
ax.hist(angles, bins=bins, density=True)
# Fiddle with labels and limits
ax.set_xlim([-np.pi/4, np.pi/4])
ax.set_xticks([-np.pi/4, 0, np.pi/4])
ax.set_xticklabels([r'$-\pi/4$', r'$0$', r'$\pi/4$'])

# bin data for our polar histogram
counts, bin = np.histogram(angles, bins=bins)
# Normalise counts to compute areas
area = counts / angles.size
# Compute corresponding radii from areas
radius = (area / np.pi)**.5

polar_ax.bar(bin[:-1], radius, align='edge', color='C0')

# Label angles according to convention
polar_ax.set_xticks([0, np.pi/4, 2*np.pi - np.pi/4])
polar_ax.set_xticklabels([r'$0$', r'$\pi/4$', r'$-\pi/4$'])

【讨论】：

伟大的贡献。我只是加快了方向统计的速度。一个开创性的参考在这里：palaeo.spb.ru/pmlibrary/pmbooks/mardia&jupp_2000.pdf.

@JayInNyc 感谢您的积极反馈 :) 您链接的文本以及 Fisher 的“循环数据的统计分析”教会了我关于方向统计的所有知识。