Mathematica - Quadratic Form Contours

Author

Quadratic Form Contours

Joerg Schaber

2005-06-03, 9:10 am

Hi,

does anybody know a simple way to calculate 2-D ellipsoids in x={x1,x2}
of a quadatric form solving xAx=b for a graphical output, i.e. contour
lines of the quadrtic form expressed as the Graphics primitive Ellipsoid?
I suppose that the option ParameterConfidenceRegion of NonlinearRegress
does something like that, but how?

best,

joerg

2005-06-07, 3:59 am

Hi Joerg,
To use "Ellipsoid" we need the semi-length and directions of the axes.
The direction we obtain from: Eigenvectors[A]
the half-length from: Sqrt[ b/Eigenvalues[A]]
Therefore, an ellipsoid would be specified by:

Ellipsoid[{x,y}, Sqrt[ b/Eigenvalues[A]], Eigenvectors[A] ]

The erason for this is:
If you turn the coordinate system in direction of the Eigenvectors, the
quadratic form looks like:

x1^2/ew1 + x2^2/ew2 == b

where ew means eigenvalue.
sincerely, Daniel

Joerg Schaber wrote:
> Hi,
>
> does anybody know a simple way to calculate 2-D ellipsoids in x={x1,x2}
> of a quadatric form solving xAx=b for a graphical output, i.e. contour
> lines of the quadrtic form expressed as the Graphics primitive Ellipsoid?
> I suppose that the option ParameterConfidenceRegion of NonlinearRegress
> does something like that, but how?
>
> best,
>
> joerg
>

Paul Abbott

2005-06-07, 9:03 am

In article <d83dqf$o1e$1@smc.vnet.net>, dh <dh@metrohm.ch> wrote:

> To use "Ellipsoid" we need the semi-length and directions of the axes.
> The direction we obtain from: Eigenvectors[A]
> the half-length from: Sqrt[ b/Eigenvalues[A]]
> Therefore, an ellipsoid would be specified by:
>
> Ellipsoid[{x,y}, Sqrt[ b/Eigenvalues[A]], Eigenvectors[A] ]

Not quite. If A is an exact matrix then Eigenvectors[A] will not be
normalised. Also, you need to transpose the Eigenvectors. The following
code does the job:

<< Statistics`MultiDescriptiveStatistics`

MatrixToEllipsoid[m_, b_, origin_:{0,0}] :=
Module[{a = Transpose[Eigenvectors[m]], r = Sqrt[b/Eigenvalues[m]]},
Ellipsoid[origin, r, a/(Norm /@ a)]]

For example, with

A = {{5, -2}, {-2, 8}}; b = 36;

the ellipsoid equation is obtained via

eqn = ExpandAll[{x, y} . A . {x, y} == b]

We can plot the contours of this as follows:

ContourPlot[First[eqn], {x, -4, 4}, {y, -4, 4},
Contours -> {b}, ContourShading -> False];

Alternatively, we can visualize the ellipse using

Show[Graphics[MatrixToEllipsoid[A,b]],
AspectRatio -> Automatic, Axes -> True]

The two pictures agree.

Cheers,
Paul
[color=darkred]
> The erason for this is:
> If you turn the coordinate system in direction of the Eigenvectors, the
> quadratic form looks like:
>
> x1^2/ew1 + x2^2/ew2 == b
>
> where ew means eigenvalue.
> sincerely, Daniel
>
> Joerg Schaber wrote:

--
Paul Abbott Phone: +61 8 6488 2734
School of Physics, M013 Fax: +61 8 6488 1014
The University of Western Australia (CRICOS Provider No 00126G)
AUSTRALIA http://physics.uwa.edu.au/~paul
http://InternationalMathematicaSymposium.org/IMS2005/

Sponsored Links