我想在 3D 图中围绕我的 PCA 集群绘制一个 2-STD 置信椭圆体； 我想要实现的是：https://www.originlab.com/fileExchange/details.aspx?fid=280 在 Python 中。我很感激任何意见。
我尝试了一些使用 Pearson 系数和集群协方差矩阵特征值的方法。 但是，椭圆体角度似乎偏离了，我无法使用以下代码重现 PC1 与 PC2 的 2D 置信椭圆：
def confidence_ellipse3D(x, y, z, ax, n_std=2.0, facecolor='none', **kwargs): """ 创建 x 和 y 的协方差置信椭圆图。
Parameters
----------
x, y : array-like, shape (n, )
Input data.
ax : matplotlib.axes.Axes
The axes object to draw the ellipse into.
n_std : float
The number of standard deviations to determine the ellipse's radiuses.
**kwargs
Forwarded to `~matplotlib.patches.Ellipse`
Returns
-------
matplotlib.patches.Ellipse
"""
if x.size != y.size or x.size != z.size:
raise ValueError("x and y must be the same size")
cov = np.cov( np.stack((x, y, z), axis=0) )
# https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5816993/#:~:text=Go%20to%3A-,PEARSON%20CORRELATION%20COEFFICIENT,%CF%81'%20as%20population%20correlation%20coefficient.
#pearsonxy = cov[0, 1]/np.sqrt(abs(cov[0, 0] * cov[1, 1])) # Formula https://slideplayer.com/10694667/37/images/slide_1.jpg
# Using a special case to obtain the eigenvalues of this
# two-dimensionl dataset.
# ell_radius_x^2 + ell_radius_y^2 = 2
_,s,rotation = linalg.svd(cov)#https://en.wikipedia.org/wiki/Covariance_matrix
radii = 1.0/np.sqrt(s)
# Calculating the stdandard deviation of x from
# the squareroot of the variance and multiplying
# with the given number of standard deviations.
scale_x = np.sqrt(abs(cov[0, 0])) * n_std
mean_x = np.mean(x)
# calculating the stdandard deviation of y ...
scale_y = np.sqrt(abs(cov[1, 1])) * n_std
mean_y = np.mean(y)
# calculating the stdandard deviation of y ...
scale_z = np.sqrt(abs(cov[2, 2])) * n_std
mean_z = np.mean(z)
# Radii corresponding to the coefficients:
radii_scaled = [radii[0]*scale_x, radii[1]*scale_y, radii[2]*scale_z]
center = [mean_x,mean_y,mean_z]
# Set of all spherical angles:
u = np.linspace(0, 2 * np.pi, 100)
v = np.linspace(0, np.pi, 100)
# Cartesian coordinates that correspond to the spherical angles:
# (this is the equation of an ellipsoid):
x = radii[0] * np.outer(np.cos(u), np.sin(v))
y = radii[1] * np.outer(np.sin(u), np.sin(v))
z = radii[2] * np.outer(np.ones_like(u), np.cos(v))
for ki in range(len(x)):
for kj in range(len(x)):
[x[ki,kj],y[ki,kj],z[ki,kj]] = np.dot([scale_x*x[ki,kj],scale_y*y[ki,kj],scale_z*z[ki,kj]], rotation) + center
# Plot:
ax.plot_surface(x, y, z, alpha=0.2)
def confidence_ellipse(x, y, ax, n_std=3.0, facecolor='none', **kwargs): """ 创建 x 和 y 的协方差置信椭圆图。
Parameters
----------
x, y : array-like, shape (n, )
Input data.
ax : matplotlib.axes.Axes
The axes object to draw the ellipse into.
n_std : float
The number of standard deviations to determine the ellipse's radiuses.
**kwargs
Forwarded to `~matplotlib.patches.Ellipse`
Returns
-------
matplotlib.patches.Ellipse
"""
if x.size != y.size:
raise ValueError("x and y must be the same size")
cov = np.cov(x, y)
pearson = cov[0, 1]/np.sqrt(abs(cov[0, 0] * cov[1, 1]))
# Using a special case to obtain the eigenvalues of this
# two-dimensionl dataset.
ell_radius_x = np.sqrt(abs(1 + pearson))
ell_radius_y = np.sqrt(abs(1 - pearson))
ellipse = Ellipse((0, 0), width=ell_radius_x * 2, height=ell_radius_y * 2,
facecolor=facecolor, **kwargs)
# Calculating the stdandard deviation of x from
# the squareroot of the variance and multiplying
# with the given number of standard deviations.
scale_x = np.sqrt(abs(cov[0, 0])) * n_std
mean_x = np.mean(x)
# calculating the stdandard deviation of y ...
scale_y = np.sqrt(abs(cov[1, 1])) * n_std
mean_y = np.mean(y)
transf = transforms.Affine2D() \
.rotate_deg(45) \
.scale(scale_x, scale_y) \
.translate(mean_x, mean_y)
ellipse.set_transform(transf + ax.transData)
return ax.add_patch(ellipse)