Skip to main content

Module ellipse

Module ellipse 

Source
Expand description

Landing ellipses: where a rocket’s landings scatter on the ground, as an ellipse that holds a chosen share of them.

Guide: Monte Carlo dispersion’s Landing ellipses section draws one from a run and says how far to trust it.

A Scatter keeps the landing points of a run (east and north of the pad, m) and the number of samples tried, so a failed flight, or one that never landed, is counted and not dropped, as in a Distribution.

§The ellipse of a normal spread

If the landings follow a two-dimensional normal distribution with mean μ and covariance Σ, the points x with (x − μ)ᵀ Σ⁻¹ (x − μ) ≤ k² fill an ellipse centered on μ. Its axes lie along the eigenvectors of Σ, its semi-axes are k √λ₁ and k √λ₂ for the eigenvalues λ₁ ≥ λ₂, and it holds the probability P(χ²₂ ≤ k²), as the left side is chi-square with two degrees of freedom. That distribution’s cumulative function is 1 − e^(−x/2), so the ellipse holding a share p (its level) has

k² = −2 ln(1 − p)

(gaussian_scale). M. Abramowitz and I. A. Stegun, Handbook of Mathematical Functions, NBS AMS 55, 1964, integrate the bivariate normal density over this ellipse to 1 − e^(−k²/2) (p. 940, eq. 26.3.21), the chi-square function with two degrees of freedom (eq. 26.4.5, p. 941). B. Wang, W. Shi and Z. Miao, “Confidence analysis of standard deviational ellipse and its extension into higher dimensional Euclidean space”, PLoS ONE 10(3), e0118537, 2015, https://doi.org/10.1371/journal.pone.0118537, derive the axes (eqs. 13–16) and the level (eqs. 19–20). For p = 50%, 90%, 95% and 99%, k² is 2 ln 2, 2 ln 10, 2 ln 20 and 4 ln 10: 1.386, 4.605, 5.991 and 9.210, as the NIST/SEMATECH e-Handbook of Statistical Methods tabulates the last three (§1.3.6.7.4, https://www.itl.nist.gov/div898/handbook/eda/section3/eda3674.htm).

For a symmetric 2 × 2 matrix Σ = [[a, b], [b, c]] (a the east variance, c the north, b their covariance) the eigenvalues are λ = (a + c)/2 ± √(((a − c)/2)² + b²) and the major axis makes the angle θ = ½ atan2(2b, a − c) with east, counter-clockwise. It is reported as a heading, clockwise from north: π/2 − θ, in [0, π). A circle (a = c, b = 0) has no major axis; its heading is reported as east’s, π/2.

Scatter::ellipse puts the sample’s mean and covariance in place of μ and Σ, its axes measured from the points themselves (Scatter::principal_axes) so that a very narrow spread keeps its width. The mean and covariance are taken on the points shifted by the first (sorted) one, the covariance with n − 1 and two passes, as Distribution’s are (T. F. Chan, G. H. Golub and R. J. LeVeque, The American Statistician 37(3), 242–247, 1983).

§The ellipse a new flight lands in

A sample’s mean and covariance are estimates, so the ellipse drawn from them holds a little less than p of the flights still to come; with few samples, much less. For normal landings the region a new flight lands in with probability exactly p, given n flights, is (x − x̄)ᵀ S⁻¹ (x − x̄) ≤ k² with

k² = 2 (n + 1)(n − 1) / (n (n − 2)) · F₂,ₙ₋₂(p) = ((n² − 1)/n) ((1 − p)^(−2/(n − 2)) − 1)

(prediction_scale, Scatter::prediction_ellipse). This follows from the new point x − x̄ being normal with covariance (1 + 1/n) Σ and independent of S, so n/(n + 1) (x − x̄)ᵀ S⁻¹ (x − x̄) is Hotelling’s T² with n − 1 degrees of freedom, which is 2(n − 1)/(n − 2) times an F with 2 and n − 2 (H. Hotelling, “The generalization of Student’s ratio”, Annals of Mathematical Statistics 2(3), 360–378, 1931, cited for the distribution and not consulted). The formula is checked against the NIST/SEMATECH e-Handbook of Statistical Methods, §6.5.4.3.4, which gives the same limit, p(m + 1)(m − 1)/(m² − mp) F(p, m − p) for p dimensions and m points, after T. P. Ryan, Statistical Methods for Quality Improvement, 2000, ch. 9 (https://www.itl.nist.gov/div898/handbook/pmc/section5/pmc5434.htm). The F distribution with 2 and m degrees of freedom has the cumulative function 1 − (1 + 2f/m)^(−m/2) (A&S eq. 26.6.4, p. 946), which inverts in closed form. As n grows, k² falls to the normal ellipse’s −2 ln(1 − p): at 200 flights and 95% it is 2.6% above it, and the semi-axes 1.3% longer.

§Whether the landings are normal

Neither ellipse is right if the landings aren’t normal, and they often aren’t: a wind whose heading is uncertain spreads them along an arc. Scatter::share_inside counts the landings an ellipse really holds. A sample that gave no landing could have landed inside or outside, so the share is a Share: a lower bound counting it outside, an upper bound counting it inside. A share far from the level means the ellipse is the wrong shape for this run; a share close to it is consistent with normal landings, not proof of them. With few landings the share runs high, as the ellipse is fitted to the same points: three points are each exactly √(4/3) standard deviations out, so even the 50% ellipse holds all three.

Every sum runs over the points sorted (east, then north), so an ellipse is bit for bit the same however the points were computed or ordered.

Structs§

Covariance
The covariance of a spread of points on the ground, m². It serializes as its three entries and reads back through Covariance::new’s checks.
Ellipse
An ellipse on the ground, holding a share of the landings. Built by Ellipse::gaussian, Scatter::ellipse or Scatter::prediction_ellipse; it serializes as its fields and reads back through checks on each.
PrincipalAxes
The axes of a Covariance: its eigenvalues, and the heading of the larger one’s eigenvector.
Scatter
Points on the ground (east and north of the pad, m), sorted, and how many samples were tried. It serializes as those two, and reads back through Scatter::new’s checks.

Functions§

gaussian_scale
The scale k of the ellipse holding the share level of a normal spread whose mean and covariance are known: k² = −2 ln(1 − level) (the module’s docs).
prediction_scale
The scale k of the ellipse a new flight lands in with probability level, from count normal landings whose mean and covariance were estimated: k² = ((n² − 1)/n) ((1 − level)^(−2/(n − 2)) − 1) (the module’s docs).