我正在尝试使用 Seaborn 在极投影上绘制双变量(联合)KDE。 Seaborn 不支持这个,Scipy 也不直接支持 angular (von Mises) KDE。
scipy gaussian_kde 和循环数据 解决了一个相关但不同的案例。相似之处是 - 随机变量是在单位圆上的线性间隔角度上定义的; KDE 被绘制。差异:我想使用 Seaborn 的 joint kernel density estimate support 来生成这种等高线图 -
但没有分类(“物种”)变化,并且在极投影上。边际地块很不错,但并不重要。
我的情况的直线版本是
import matplotlib
import pandas as pd
from matplotlib import pyplot as plt
import numpy as np
import seaborn as sns
from numpy.random._generator import default_rng
angle = np.repeat(
np.deg2rad(
np.arange(0, 360, 10)
),
100,
)
rand = default_rng(seed=0)
data = pd.Series(
rand.normal(loc=50, scale=10, size=angle.size),
index=pd.Index(angle, name='angle'),
name='power',
)
matplotlib.use(backend='TkAgg')
joint = sns.JointGrid(
data.reset_index(),
x='angle', y='power'
)
joint.plot_joint(sns.kdeplot, bw_adjust=0.7, linewidths=1)
plt.show()
但这是在错误的投影中显示的,并且0和360之间的角度也不应该有递减的轮廓线。
当然,正如 Creating a circular density plot using matplotlib and seaborn 解释的那样,在极投影中使用现有高斯 KDE 的天真方法是无效的,即使我想我也做不到,因为
axisgrid.py
在没有参数的情况下对子图设置进行硬编码:
f = plt.figure(figsize=(height, height))
gs = plt.GridSpec(ratio + 1, ratio + 1)
ax_joint = f.add_subplot(gs[1:, :-1])
ax_marg_x = f.add_subplot(gs[0, :-1], sharex=ax_joint)
ax_marg_y = f.add_subplot(gs[1:, -1], sharey=ax_joint)
我开始使用猴子修补方法:
import scipy.stats._kde
import numpy as np
def von_mises_estimate(
points: np.ndarray,
values: np.ndarray,
xi: np.ndarray,
cho_cov: np.ndarray,
dtype: np.dtype,
real: int = 0
) -> np.ndarray:
"""
Mimics the signature of gaussian_kernel_estimate
https://github.com/scipy/scipy/blob/main/scipy/stats/_stats.pyx#L740
"""
# https://stackoverflow.com/a/44783738
# Will make this a parameter
kappa = 20
# I am unclear on how 'values' would be used here
class VonMisesKDE(scipy.stats._kde.gaussian_kde):
def __call__(self, points: np.ndarray) -> np.ndarray:
points = np.atleast_2d(points)
result = von_mises_estimate(
self.dataset.T,
self.weights[:, None],
points.T,
self.inv_cov,
points.dtype,
)
return result[:, 0]
import seaborn._statistics
seaborn._statistics.gaussian_kde = VonMisesKDE
这确实成功地调用了默认的高斯函数,但是(1)它是不完整的,并且(2)我不清楚是否有可能说服联合绘图方法使用新的投影。
一个非常扭曲和低质量的预览,通过 Gimp 转换:
虽然径向轴会增加而不是从中心向外减少。
我过去通过结合使用
seaborn
和matplotlib
的极投影来做到这一点。这是一个例子:
import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns
n_data = 1000
data_phase = np.random.rand(n_data) * 1.2 * np.pi
data_amp = np.random.randn(n_data)
fig = plt.figure()
ax = fig.add_subplot(111, projection='polar')
ax.scatter(data_phase, data_amp, vmin=0, vmax=2 * np.pi, s=10, alpha=0.3)
ax.set_thetagrids(angles=np.linspace(0, 360, 5));
sns.kdeplot(x=data_phase, y=data_amp, n_levels=5, c='k', ax=ax)
希望您能够从那里调整它以满足您的需求?
这里有一个方法的想法:
sin
和cos
jointplot
(或kdeplot
,这可以包括hue
)from matplotlib import pyplot as plt
import seaborn as sns
import pandas as pd
import numpy as np
# test data from https://www.kaggle.com/datasets/muthuj7/weather-dataset
df = pd.read_csv('weatherhistory.csv')[['Wind Speed (km/h)', 'Wind Bearing (degrees)']].rename(
columns={'Wind Bearing (degrees)': 'angle', 'Wind Speed (km/h)': 'power'})
df['angle'] = np.radians(df['angle'])
df['x'] = df['power'] * np.cos(df['angle'])
df['y'] = df['power'] * np.sin(df['angle'])
fig = plt.figure(figsize=(10, 10))
grid_ratio = 5
gs = plt.GridSpec(grid_ratio + 1, grid_ratio + 1)
ax_joint = fig.add_subplot(gs[1:, :-1])
ax_marg_x = fig.add_subplot(gs[0, :-1])
ax_marg_y = fig.add_subplot(gs[1:, -1])
sns.kdeplot(data=df, x='x', y='y', bw_adjust=0.7, linewidths=1, ax=ax_joint)
ax_joint.set_aspect('equal', adjustable='box') # equal aspect ratio is needed for a polar plot
ax_joint.axis('off')
xmin, xmax = ax_joint.get_xlim()
xrange = max(-xmin, xmax)
ax_joint.set_xlim(-xrange, xrange) # force 0 at center
ymin, ymax = ax_joint.get_ylim()
yrange = max(-ymin, ymax)
ax_joint.set_ylim(-yrange, yrange) # force 0 at center
ax_polar = fig.add_subplot(projection='polar')
ax_polar.set_facecolor('none') # make transparent
ax_polar.set_position(pos=ax_joint.get_position())
ax_polar.set_rlim(0, max(xrange, yrange))
# add kdeplot of power as marginal y
sns.kdeplot(y=df['power'], ax=ax_marg_y)
ax_marg_y.set_ylim(0, df['power'].max() * 1.1)
ax_marg_y.set_xlabel('')
ax_marg_y.set_ylabel('')
ax_marg_y.text(1, 0.5, 'power', transform=ax_marg_y.transAxes, ha='center', va='center')
sns.despine(ax=ax_marg_y, bottom=True)
# add kdeplot of angles as marginal x, repeat the angles shifted -360 and 360 degrees to enable wrap-around
angles = np.degrees(df['angle'])
angles_trippled = np.concatenate([angles - 360, angles, angles + 360])
sns.kdeplot(x=angles_trippled, ax=ax_marg_x)
ax_marg_x.set_xlim(0, 360)
ax_marg_x.set_xticks(np.arange(0, 361, 45))
ax_marg_x.set_xlabel('')
ax_marg_x.set_ylabel('')
ax_marg_x.text(0.5, 1, 'angle', transform=ax_marg_x.transAxes, ha='center', va='center')
sns.despine(ax=ax_marg_x, left=True)
plt.show()
PS:这就是填充版本的样子(带有
cmap='turbo'
):
如果你想要0在顶部,让角度顺时针旋转,你需要在调用2D
x=
时切换y=
和kdeplot
。
sns.kdeplot(data=df, x='y', y='x', bw_adjust=0.7, fill=True, cmap='turbo', ax=ax_joint)
# ....
ax_polar = fig.add_subplot(projection='polar')
ax_polar.set_theta_zero_location('N')
ax_polar.set_theta_direction('clockwise')