blocks|key|189048|text|您可以尝试对数据矩阵进行奇异值分解。https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.linalg.svd.html|type|unstyled|depth|inlineStyleRanges|entityRanges|offset|length|data|189049|首先，从每列中减去X，Y的平均值，对数据进行中心化处理。|189050|X=X-np.mean(X)
Y=Y-np.mean(Y)

D=np.vstack(X,Y)|code-block|syntax|javascript|189051|然后应用奇异值分解，提取-eigenvalues+(s的成员)+->轴长-eigenvectors(U)+->轴方向。|189052|U,+s,+V+=+np.linalg.svd(D,+full_matrices=True)|189053|这应该是最不合适的。当然，事情可能会变得更复杂，请参阅https://www.emis.de/journals/BBMS/Bulletin/sup962/gander.pdf|189054|entityMap|0|LINK|mutability|MUTABLE|url|https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.linalg.svd.html|1|https://www.emis.de/journals/BBMS/Bulletin/sup962/gander.pdf^0|I|29|0|0|0|0|0|0|R|1O|1|0^^$0|@$1|2|3|4|5|6|7|10|8|@]|9|@$A|11|B|12|1|13]]|C|$]]|$1|D|3|E|5|6|7|14|8|@]|9|@]|C|$]]|$1|F|3|G|5|H|7|15|8|@]|9|@]|C|$I|J]]|$1|K|3|L|5|6|7|16|8|@]|9|@]|C|$]]|$1|M|3|N|5|H|7|17|8|@]|9|@]|C|$I|J]]|$1|O|3|P|5|6|7|18|8|@]|9|@$A|19|B|1A|1|1B]]|C|$]]|$1|Q|3|-4|5|6|7|1C|8|@]|9|@]|C|$]]]|R|$S|$5|T|U|V|C|$W|X]]|Y|$5|T|U|V|C|$W|Z]]]]

You could try a Singular Value Decomposition of the data matrix.
<a href="https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.linalg.svd.html" rel="nofollow noreferrer">https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.linalg.svd.html</a>

First center the data by subtracting mean values of X,Y from each column respectively.

<pre><code>X=X-np.mean(X)
Y=Y-np.mean(Y)

D=np.vstack(X,Y)
</code></pre>

Then, apply SVD and extract
-eigenvalues (members of s) -> axis length
-eigenvectors(U) -> axis orientation

<pre><code>U, s, V = np.linalg.svd(D, full_matrices=True)
</code></pre>

This should be a least-squares fit.
Of course, things can get more complicated than this, please see
<a href="https://www.emis.de/journals/BBMS/Bulletin/sup962/gander.pdf" rel="nofollow noreferrer">https://www.emis.de/journals/BBMS/Bulletin/sup962/gander.pdf</a>

blocks|key|3640282|text|按照ErroriSalvo的建议，下面是用奇异值分解拟合椭圆的完整过程。数组x，y是给定点的坐标，假设有N个点。然后由形状(2，N)的中心坐标阵的SVD得到U，S，V。因此，U是2×2正交矩阵(旋转)，S是长度为2的向量(奇异值)，V是N乘N的正交矩阵。|type|unstyled|depth|inlineStyleRanges|entityRanges|data|3640283|将单位圆转化为最佳拟合椭圆的线性映射是|3640284|sqrt(2/N)+*+U+*+diag(S)+|code-block|syntax|javascript|3640285|其中diag(S)是对角线上具有奇异值的对角线矩阵。为了了解为什么需要sqrt(2/N)的因子，假设点x，y是从单位圆中均匀取出来的。然后sum(x**2)+%2B+sum(y**2)是N，因此坐标矩阵由长度为sqrt(N/2)的两个正交行组成，因此它的范数(最大奇异值)是sqrt(N/2)。我们需要把这个降到1，才能得到单位圆。|offset|length|style|CODE|3640286|N+=+300
t+=+np.linspace(0,+2*np.pi,+N)
x+=+5*np.cos(t)+%2B+0.2*np.random.normal(size=N)+%2B+1
y+=+4*np.sin(t%2B0.5)+%2B+0.2*np.random.normal(size=N)
plt.plot(x,+y,+'.')+++++#+given+points

xmean,+ymean+=+x.mean(),+y.mean()
x+-=+xmean
y+-=+ymean
U,+S,+V+=+np.linalg.svd(np.stack((x,+y)))

tt+=+np.linspace(0,+2*np.pi,+1000)
circle+=+np.stack((np.cos(tt),+np.sin(tt)))++++#+unit+circle
transform+=+np.sqrt(2/N)+*+U.dot(np.diag(S))+++#+transformation+matrix
fit+=+transform.dot(circle)+%2B+np.array([[xmean],+[ymean]])
plt.plot(fit[0,+:],+fit[1,+:],+'r')
plt.show()|3640287|​|3640288|📷|atomic|3640289|3640290|但是如果你假设没有旋转，那么np.sqrt(2/N)+*+S就是你所需要的，这就是椭圆方程中的a和b。|3640291|entityMap|0|IMAGE|mutability|IMMUTABLE|imageUrl|https://developer.qcloudimg.com/http-save/yehe-900000/5a1933ceb256d165b0c78ee37bcba4b9.png|imageAlt^0|0|0|0|1X|L|0|0|0|0|1|0|0|0|E|G|1B|1|1D|1|0^^$0|@$1|2|3|4|5|6|7|17|8|@]|9|@]|A|$]]|$1|B|3|C|5|6|7|18|8|@]|9|@]|A|$]]|$1|D|3|E|5|F|7|19|8|@]|9|@]|A|$G|H]]|$1|I|3|J|5|6|7|1A|8|@$K|1B|L|1C|M|N]]|9|@]|A|$]]|$1|O|3|P|5|F|7|1D|8|@]|9|@]|A|$G|H]]|$1|Q|3|R|5|6|7|1E|8|@]|9|@]|A|$]]|$1|S|3|T|5|U|7|1F|8|@]|9|@$K|1G|L|1H|1|1I]]|A|$]]|$1|V|3|R|5|6|7|1J|8|@]|9|@]|A|$]]|$1|W|3|X|5|6|7|1K|8|@$K|1L|L|1M|M|N]|$K|1N|L|1O|M|N]|$K|1P|L|1Q|M|N]]|9|@]|A|$]]|$1|Y|3|-4|5|6|7|1R|8|@]|9|@]|A|$]]]|Z|$10|$5|11|12|13|A|$14|15|16|-4]]]]

Following the suggestion by ErroriSalvo, here is the complete process of fitting an ellipse using the SVD. The arrays x, y are coordinates of the given points, let's say there are N points. Then U, S, V are obtained from the SVD of the centered coordinate array of shape (2, N). So, U is a 2 by 2 orthogonal matrix (rotation), S is a vector of length 2 (singular values), and V, which we do not need, is an N by N orthogonal matrix. 

The linear map transforming the unit circle to the ellipse of best fit is 

<pre><code>sqrt(2/N) * U * diag(S) 
</code></pre>

where diag(S) is the diagonal matrix with singular values on the diagonal. To see why the factor of sqrt(2/N) is needed, imagine that the points x, y are taken uniformly from the unit circle. Then <code>sum(x**2) + sum(y**2)</code> is N, and so the coordinate matrix consists of two orthogonal rows of length sqrt(N/2), hence its norm (the largest singular value) is sqrt(N/2). We need to bring this down to 1 to have the unit circle.

<pre><code>N = 300
t = np.linspace(0, 2*np.pi, N)
x = 5*np.cos(t) + 0.2*np.random.normal(size=N) + 1
y = 4*np.sin(t+0.5) + 0.2*np.random.normal(size=N)
plt.plot(x, y, '.') # given points

xmean, ymean = x.mean(), y.mean()
x -= xmean
y -= ymean
U, S, V = np.linalg.svd(np.stack((x, y)))

tt = np.linspace(0, 2*np.pi, 1000)
circle = np.stack((np.cos(tt), np.sin(tt))) # unit circle
transform = np.sqrt(2/N) * U.dot(np.diag(S)) # transformation matrix
fit = transform.dot(circle) + np.array([[xmean], [ymean]])
plt.plot(fit[0, :], fit[1, :], 'r')
plt.show()
</code></pre>

<a href="https://i.stack.imgur.com/p7c9n.png" rel="noreferrer"><img src="https://i.stack.imgur.com/p7c9n.png" alt="example"></a>

But if you assume that there is no rotation, then <code>np.sqrt(2/N) * S</code> is all you need; these are <code>a</code> and <code>b</code> in the equation of the ellipse.

blocks|key|3640342|text|这可以直接用最小二乘来解决。你可以用最小化数量平方之和(α*+x_i%5E2+%2Bβ*+y_i%5E2+-+1)来表示，其中α为1/a%5E2，β为1/b%5E2。所有的x_i在X中，y_i在Y中，因此可以找到x_i的极小值，其中A是Nx2矩阵(即X%5E2，Y%5E2)，x是列向量α；β和b是所有矩阵的列向量。|type|unstyled|depth|inlineStyleRanges|entityRanges|data|3640343|下面的代码解决了形式为Ax%5E2+%2B+Bxy+%2B+Cy%5E2+%2B+Dx+%2BEy+=1的椭圆的更一般问题，尽管思路完全相同。打印语句给出0.0776x%5E2+%2B+0.0315xy%2B0.125y%5E2%2B0.00457x%2B0.00314y+=1，生成的椭圆的图像也在下面|3640344|import+numpy+as+np
import+matplotlib.pyplot+as+plt
alpha+=+5
beta+=+3
N+=+500
DIM+=+2

np.random.seed(2)

#+Generate+random+points+on+the+unit+circle+by+sampling+uniform+angles
theta+=+np.random.uniform(0,+2*np.pi,+(N,1))
eps_noise+=+0.2+*+np.random.normal(size=[N,1])
circle+=+np.hstack([np.cos(theta),+np.sin(theta)])

#+Stretch+and+rotate+circle+to+an+ellipse+with+random+linear+tranformation
B+=+np.random.randint(-3,+3,+(DIM,+DIM))
noisy_ellipse+=+circle.dot(B)+%2B+eps_noise

#+Extract+x+coords+and+y+coords+of+the+ellipse+as+column+vectors
X+=+noisy_ellipse[:,0:1]
Y+=+noisy_ellipse[:,1:]

#+Formulate+and+solve+the+least+squares+problem+%7C%7CAx+-+b+%7C%7C%5E2
A+=+np.hstack([X**2,+X+*+Y,+Y**2,+X,+Y])
b+=+np.ones_like(X)
x+=+np.linalg.lstsq(A,+b)[0].squeeze()

#+Print+the+equation+of+the+ellipse+in+standard+form
print('The+ellipse+is+given+by+{0:.3}x%5E2+%2B+{1:.3}xy%2B{2:.3}y%5E2%2B{3:.3}x%2B{4:.3}y+=+1'.format(x[0],+x[1],x[2],x[3],x[4]))

#+Plot+the+noisy+data
plt.scatter(X,+Y,+label='Data+Points')

#+Plot+the+original+ellipse+from+which+the+data+was+generated
phi+=+np.linspace(0,+2*np.pi,+1000).reshape((1000,1))
c+=+np.hstack([np.cos(phi),+np.sin(phi)])
ground_truth_ellipse+=+c.dot(B)
plt.plot(ground_truth_ellipse[:,0],+ground_truth_ellipse[:,1],+'k--',+label='Generating+Ellipse')

#+Plot+the+least+squares+ellipse
x_coord+=+np.linspace(-5,5,300)
y_coord+=+np.linspace(-5,5,300)
X_coord,+Y_coord+=+np.meshgrid(x_coord,+y_coord)
Z_coord+=+x[0]+*+X_coord+**+2+%2B+x[1]+*+X_coord+*+Y_coord+%2B+x[2]+*+Y_coord**2+%2B+x[3]+*+X_coord+%2B+x[4]+*+Y_coord
plt.contour(X_coord,+Y_coord,+Z_coord,+levels=[1],+colors=('r'),+linewidths=2)

plt.legend()
plt.xlabel('X')
plt.ylabel('Y')
plt.show()|code-block|syntax|javascript|3640345|​|3640346|📷|atomic|offset|length|3640347|3640348|entityMap|0|IMAGE|mutability|IMMUTABLE|imageUrl|https://developer.qcloudimg.com/http-save/yehe-900000/539bd04c58805ae9c3e46b07dcf667b5.png|imageAlt^0|0|0|0|0|0|1|0|0|0^^$0|@$1|2|3|4|5|6|7|Z|8|@]|9|@]|A|$]]|$1|B|3|C|5|6|7|10|8|@]|9|@]|A|$]]|$1|D|3|E|5|F|7|11|8|@]|9|@]|A|$G|H]]|$1|I|3|J|5|6|7|12|8|@]|9|@]|A|$]]|$1|K|3|L|5|M|7|13|8|@]|9|@$N|14|O|15|1|16]]|A|$]]|$1|P|3|J|5|6|7|17|8|@]|9|@]|A|$]]|$1|Q|3|-4|5|6|7|18|8|@]|9|@]|A|$]]]|R|$S|$5|T|U|V|A|$W|X|Y|-4]]]]

This can be solved directly using least squares. You can frame this as minimizing the sum of squares of quantity (alpha * x_i^2 + beta * y_i^2 - 1) where alpha is 1/a^2 and beta is 1/b^2. You have all the x_i's in X and the y_i's in Y so you can find the minimizer of ||Ax - b||^2 where A is an Nx2 matrix (i.e. [X^2, Y^2]), x is the column vector [alpha; beta] and b is column vector of all ones.

The following code solves the more general problem for an ellipse of the form Ax^2 + Bxy + Cy^2 + Dx +Ey = 1 though the idea is exactly the same. The print statement gives 0.0776x^2 + 0.0315xy+0.125y^2+0.00457x+0.00314y = 1 and the image of the ellipse generated is also below

<pre><code>import numpy as np
import matplotlib.pyplot as plt
alpha = 5
beta = 3
N = 500
DIM = 2

np.random.seed(2)

# Generate random points on the unit circle by sampling uniform angles
theta = np.random.uniform(0, 2*np.pi, (N,1))
eps_noise = 0.2 * np.random.normal(size=[N,1])
circle = np.hstack([np.cos(theta), np.sin(theta)])

# Stretch and rotate circle to an ellipse with random linear tranformation
B = np.random.randint(-3, 3, (DIM, DIM))
noisy_ellipse = circle.dot(B) + eps_noise

# Extract x coords and y coords of the ellipse as column vectors
X = noisy_ellipse[:,0:1]
Y = noisy_ellipse[:,1:]

# Formulate and solve the least squares problem ||Ax - b ||^2
A = np.hstack([X**2, X * Y, Y**2, X, Y])
b = np.ones_like(X)
x = np.linalg.lstsq(A, b)[0].squeeze()

# Print the equation of the ellipse in standard form
print('The ellipse is given by {0:.3}x^2 + {1:.3}xy+{2:.3}y^2+{3:.3}x+{4:.3}y = 1'.format(x[0], x[1],x[2],x[3],x[4]))

# Plot the noisy data
plt.scatter(X, Y, label='Data Points')

# Plot the original ellipse from which the data was generated
phi = np.linspace(0, 2*np.pi, 1000).reshape((1000,1))
c = np.hstack([np.cos(phi), np.sin(phi)])
ground_truth_ellipse = c.dot(B)
plt.plot(ground_truth_ellipse[:,0], ground_truth_ellipse[:,1], 'k--', label='Generating Ellipse')

# Plot the least squares ellipse
x_coord = np.linspace(-5,5,300)
y_coord = np.linspace(-5,5,300)
X_coord, Y_coord = np.meshgrid(x_coord, y_coord)
Z_coord = x[0] * X_coord ** 2 + x[1] * X_coord * Y_coord + x[2] * Y_coord**2 + x[3] * X_coord + x[4] * Y_coord
plt.contour(X_coord, Y_coord, Z_coord, levels=[1], colors=('r'), linewidths=2)

plt.legend()
plt.xlabel('X')
plt.ylabel('Y')
plt.show()
</code></pre>

<a href="https://i.stack.imgur.com/xTZTh.png" rel="noreferrer"><img src="https://i.stack.imgur.com/xTZTh.png" alt="enter image description here"></a>

I would like to fit a 2D array by an elliptic function: (x / a)² + (y / b)² = 1 ----> (and so get the a and b)

And then, be able to replot it on my graph. 
I found many examples on internet, but no one with this simple Cartesian equation. I probably have searched badly ! I think a basic solution for this problem could help many people. 

Here is an example of the data: 

<a href="https://i.stack.imgur.com/ytPzP.png" rel="noreferrer"><img src="https://i.stack.imgur.com/ytPzP.png" alt="enter image description here"></a>

Sadly, I can not put the values... So let's assume that I have an X,Y arrays defining the coordinates of each of those points.

How to fit a 2D ellipse to given points

翻译质量差，导致语言生硬或混乱。

没有提供实际的解决方法或示例。

解答不清晰，无法理解或解决问题。

页面排版不美观，阅读体验差。

文章

问答

视频

教程

学习中心

腾讯云实验室

直播

竞赛

腾讯云代码分析专区

腾讯iOA零信任安全管理系统专区

腾讯云架构师技术同盟交流圈

腾讯云数据库专区

腾讯云智能顾问专区

腾讯云原生专区

腾讯混元专区

腾讯云TCE专区

腾讯云Lighthouse专区

腾讯云HAI专区

腾讯云Edgeone专区

腾讯云存储专区

腾讯云智能专区

腾讯轻联专区 

腾讯云开发专区

TAPD专区

腾讯轻量云游戏服专区

EdgeOne AI 安全实战专区

腾讯云最具价值专家

腾讯云架构师技术同盟

腾讯云创作之星

腾讯云开发者先锋

腾讯云代码助手

云原生构建

TAPD 敏捷项目管理

Cloud Studio

SDK中心

API中心

命令行工具

涵盖代码开发、场景应用、自动测试全流程，助你从零构建专属AI助手

一站式MCP教程库，解锁AI应用新玩法

聚焦“写作效率、视觉美观与运行性能”三方面进行全面升级，为您提供更高效、稳定的创作环境

社区富文本&Markdown编辑器全新改版上线，欢迎大家体验!

诚挚邀请您参与本次调研，分享您的真实使用感受与建议。您的反馈至关重要，感谢您的支持与参与！

社区新版编辑器体验调研

我想用椭圆函数来拟合二维阵列：(x /a)2+ (y /b)2= 1 ->(因此得到a和b)。然后，可以在我的图表上重绘它。我在网上找到了很多例子，但是没有人用这个简单的笛卡尔方程。我可能搜索得很糟糕！我认为解决这个问题的基本办法可以帮助很多人。以下是数据的一个示例：​​可悲的是，我不能把这些价值观..。让我们假设，我...

问二维椭圆与给定点的拟合方法
EN

Stack Overflow用户

社区

活动

圈层

关于

腾讯云开发者

热门产品

热门推荐

更多推荐

问二维椭圆与给定点的拟合方法EN