9
$\begingroup$

My problem is as follows: I drop 40 balls at once from a certain point, a few meters over the floor. The balls roll, and comes to a rest. Using computer vision, I calculate the center of mass in the X-Y plane. I am only interested in the distance from the center of mass to each ball, which is calculated using simple geometry. Now, I want to know the one-sided standard deviation from the center. So, I would be able to know that a certain number of balls are within one std radius, more balls within 2*std radius and so on. How do I calculate the one-sided standard deviation? A normal approach would state that half of the balls be on the "negative side" of 0 mean. This of course makes no sense in this experiment. Do I have to make sure that the balls conform to the standard distribution? Thank you for any help.

$\endgroup$
13
$\begingroup$

To characterize the amount of 2D dispersion around the centroid, you just want the (root) mean squared distance,

$$\hat\sigma=\text{RMS} = \sqrt{\frac{1}{n}\sum_i\left((x_i - \bar{x})^2 + (y_i - \bar{y})^2\right)}.$$

In this formula, $(x_i, y_i), i=1, 2, \ldots, n$ are the point coordinates and their centroid (point of averages) is $(\bar{x}, \bar{y}).$


The question asks for the distribution of the distances. When the balls have an isotropic bivariate Normal distribution around their centroid--which is a standard and physically reasonable assumption--the squared distance is proportional to a chi-squared distribution with two degrees of freedom (one for each coordinate). This is a direct consequence of one definition of the chi-squared distribution as a sum of squares of independent standard normal variables, because $$x_i - \bar{x} = \frac{n-1}{n}x_i - \sum_{j\ne i}\frac{1}{n}x_j$$ is a linear combination of independent normal variates with expectation $$\mathbb{E}[x_i - \bar{x}] = \frac{n-1}{n}\mathbb{E}[x_i] -\sum_{j\ne i}\frac{1}{n}\mathbb{E}[x_j] = 0.$$ Writing the common variance of the $x_i$ as $\sigma^2$, $$\mathbb{E}[\left(x_i -\bar{x}\right)^2]=\text{Var}(x_i - \bar{x}) = \left(\frac{n-1}{n}\right)^2\text{Var}(x_i) + \sum_{j\ne i}\left(\frac{1}{n}\right)^2\text{Var}(x_j) = \frac{n-1}{n}\sigma^2.$$ The assumption of anisotropy is that the $y_j$ have the same distribution as the $x_i$ and are independent of them, so an identical result holds for the distribution of $(y_j - \bar{y})^2$. This establishes the constant of proportionality: the squares of the distances have a chi-squared distribution with two degrees of freedom, scaled by $\frac{n-1}{n}\sigma^2$.

The most severe test of these equations is the case $n=2$, for then the fraction $\frac{n-1}{n}$ differs the most from $1$. By simulating the experiment, both for $n=2$ and $n=40$, and overplotting the histograms of squared distances with the scaled chi-squared distributions (in red), we can verify this theory.

Figure

Each row shows the same data: on the left the x-axis is logarithmic; on the right it shows the actual squared distance. The true value of $\sigma$ for these simulations was set to $1$.

These results are for 100,000 iterations with $n=2$ and 50,000 iterations with $n=40$. The agreements between the histograms and chi-squared densities are excellent.


Although $\sigma^2$ is unknown, it can be estimated in various ways. For instance, the mean squared distance should be $\frac{n-1}{n}\sigma^2$ times the mean of $\chi^2_2$, which is $2$. With $n=40$, for example, estimate $\sigma^2$ as $\frac{40}{39}/2$ times the mean squared distance. Thus an estimate of $\sigma$ would be $\sqrt{40/78}$ times the RMS distance. Using values of the $\chi^2_2$ distribution we can then say that:

  • Approximately 39% of the distances will be less than $\sqrt{39/40}\hat\sigma$, because 39% of a $\chi^2_2$ distribution is less than $1$.

  • Approximately 78% of the distances will be less than $\sqrt{3}$ times $\sqrt{39/40}\hat\sigma$, because 78% of a $\chi^2_2$ distribution is less than $3$.

And so on, for any multiple you care to use in place of $1$ or $3$. As a check, in the simulations for $n=40$ plotted previously, the actual proportions of squared distances less than $1, 2, \ldots, 10$ times $\frac{n-1}{n}\hat\sigma^2$ were

0.3932 0.6320 0.7767 0.8647 0.9178 0.9504 0.9700 0.9818 0.9890 0.9933

The theoretical proportions are

0.3935 0.6321 0.7769 0.8647 0.9179 0.9502 0.9698 0.9817 0.9889 0.9933

The agreement is excellent.


Here is R code to conduct and analyze the simulations.

f <- function(n, n.iter, x.min=0, x.max=Inf, plot=TRUE) {
  #
  # Generate `n.iter` experiments in which `n` locations are generated using
  # standard normal variates for their coordinates.
  #
  xy <- array(rnorm(n*2*n.iter), c(n.iter,2,n))
  #
  # Compute the squared distances to the centers for each experiment.
  #
  xy.center <- apply(xy, c(1,2), mean)
  xy.distances2 <- apply(xy-array(xy.center, c(n.iter,2,n)), c(1,3), 
                         function(z) sum(z^2))
  #
  # Optionally plot histograms.
  #
  if(plot) {
    xy.plot <- xy.distances2[xy.distances2 >= x.min & xy.distances2 <= x.max]

    hist(log(xy.plot), prob=TRUE, breaks=30,
         main=paste("Histogram of log squared distance, n=", n),
         xlab="Log squared distance")
    curve(dchisq(n/(n-1) * exp(x), df=2) * exp(x) * n/(n-1), 
          from=log(min(xy.plot)), to=log(max(xy.plot)), 
          n=513, add=TRUE, col="Red", lwd=2)

    hist(xy.plot, prob=TRUE, breaks=30,
         main=paste("Histogram of squared distance, n=", n),
         xlab="Squared distance")
    curve(n/(n-1) * dchisq(n/(n-1) * x, df=2), 
          from=min(xy.plot), to=max(xy.plot), 
          n=513, add=TRUE, col="Red", lwd=2)  
  }
  return(xy.distances2)
}
#
# Plot the histograms and compare to scaled chi-squared distributions.
#
par(mfrow=c(2,2))
set.seed(17)
xy.distances2 <- f(2, 10^5, exp(-6), 6)
xy.distances2 <- f(n <- 40, n.iter <- 50000, exp(-6), 12)
#
# Compare the last simulation to cumulative chi-squared distributions.
#
sigma.hat <- sqrt((n / (2*(n-1)) * mean(xy.distances2)))
print(cumsum(tabulate(cut(xy.distances2, 
                    (0:10) * (n-1)/n * sigma.hat^2))) / (n*n.iter), digits=4)
print(pchisq(1:10, df=2), digits=4)
$\endgroup$
  • 2
    $\begingroup$ Thank you for a very comprehensive answer. I cannot quite understand how the RMS formula can describe standard deviation without dividing by the number of balls. If you compare it to http://en.wikipedia.org/wiki/Root-mean-square_deviation_(bioinformatics they have divided the sum by N. Should the sum be divided by N or N-1 (since 40 balls is just a selection from a population of balls?) $\endgroup$ – K_scheduler Aug 1 '13 at 9:57
  • $\begingroup$ After doing the calculations again, it seems like sqrt(SDx^2+SDy^2) is what I am after. This will give me a radius for a circle that contains all balls with a probability of 65%, right? $\endgroup$ – K_scheduler Aug 1 '13 at 12:10
  • $\begingroup$ That's an equivalent formula for the RMS, but the 65% value is incorrect, as explained in this answer. $\endgroup$ – whuber Aug 1 '13 at 12:20
  • 2
    $\begingroup$ @nali All those points are clearly made in my answer here. $\endgroup$ – whuber Sep 13 '15 at 13:12
  • 4
    $\begingroup$ @nali Your posts here go beyond the bounds of propriety in their rudeness and ad hominem attacks. Although I am not worried about being considered ignorant or stupid, as a moderator of this site I have to be concerned about keeping the discourse civil and therefore cannot tolerate the vituperation you are posting. Accordingly, I have deleted your latest comment. If I see comments from you that are similarly rude, toward anyone at all, I will delete them without further notice and I (or other moderators) will take immediate steps to limit your interactions on this site. $\endgroup$ – whuber Sep 13 '15 at 17:09
4
$\begingroup$

I think you have some things a bit confused. It's true that distance can't be negative, but that doesn't affect calculation of the standard deviation. Although it means the distribution of distances can't be exactly normal, it could still be close; but even if its is far from normal, there is still a standard deviation.

Also, there is no "one sided" standard deviation - you may be thinking of hypothesis tests (which can be one sided or two sided). In your title, you say mean is 0, but the mean distance won't be 0 (unless the balls are in a stack 40 balls high!) and you say there are limits - there could be limits, if the balls are dropped in a room then they can't be farther from the center than the distance to the nearest wall. But unless some of the balls bounce against a wall, that won't affect things.

So, once you have the 40 distances you calculate the standard deviation (and mean, median, interquartile range, etc) using standard methods. You can also make plots of the distance (e.g. quantile normal plot, box plot) to see if it is roughly normally distributed (if that's of interest).

$\endgroup$
  • $\begingroup$ Thank you Peter, I did not express myself correctly. Let me try to clarify: Imagine the scene from above. You calculate the mean distance, it will be illustrated as a circle around the center of mass (mean distance = radius). Now, +/- std deviation from this will yield a smaller circle and a larger circle. I don't want to know the standard deviation of the mean distance to center of mass, but rather the standard deviation from the center of mass outward. In other words, within what radius from the center of mass is 68.2% (one standard deviation) of the balls situated. $\endgroup$ – K_scheduler Jul 26 '13 at 12:09
  • $\begingroup$ Oh, OK. Then I think this is not a statistics problem but a math problem; finding where 68.2% will fall is known... I forget the answer but it involve $\pi$. $\endgroup$ – Peter Flom Jul 26 '13 at 12:53
  • $\begingroup$ You might be right in your first answer. From what i have found, using the radial standard deviation should do the trick. RSD = sqrt(SDx^2 + SDy^2) $\endgroup$ – K_scheduler Jul 26 '13 at 13:04
0
$\begingroup$

Its been a while since this was asked, but the answer to the question is that this is the 2D distribution named the Rayleigh distribution. Here the assumption is that the Rayleigh shape factor is equal to both the standard deviations of the X and Y coordinates. In practice the value of the shape factor would be calculated from the pooled average of the standard deviation of X and Y.

starting with $$ X \sim \mathcal{N}(\mu_x,\sigma_x^2)$$, and $$Y \sim \mathcal{N}(\mu_y,\sigma_y^2)$$

use bivariant normal distribution. $$ f(x,y) = \frac{1}{2 \pi \sigma_x \sigma_y \sqrt{1-\rho^2}} \exp\left( -\frac{1}{2(1-\rho^2)}\left[ \frac{(x-\mu_x)^2}{\sigma_x^2} + \frac{(y-\mu_y)^2}{\sigma_y^2} - \frac{2\rho(x-\mu_x)(y-\mu_y)}{\sigma_x \sigma_y} \right] \right)$$

translate to point $$(\mu_x, \mu_y)$$ and assume $$\rho = 0$$.

Also assume that $$\sigma_x^2 = \sigma_y^2$$ so replace both with $$\sigma^2$$

then the 2-D distribution is expressed as the radius around point $$(\mu_x, \mu_y)$$ which is known as the Rayleigh distribution.

$$PDF(r; \sigma) = \frac{r}{\sigma^2 } \exp\left( - \frac{r^2}{2\sigma^2} \right) $$ where $$\sigma = \sigma_x = \sigma_y$$ and $$r_i = \sqrt{(x_i - \mu_x)^2 + (y_i - \mu_y)^2}$$

$$ CDF(r; \sigma) = 1 - \exp\left( - \frac{r^2}{2\sigma^2} \right)$$

Of course this is for the continuous distribution. For a sample of just 40 balls there is no exact solution. You'd need to do a Monte Carlo Analysis with a sample of 40 balls. Taylor, M. S. & Grubbs, Frank E. (1975). "Approximate Probability Distributions for the Extreme Spread" found estimates for the Chi distribution and the log-normal for that would fit the distribution of a sample.


Edit - Despite Wuber's doubt, the theoretical proportions he calculated are:

0.3935 0.6321 0.7769 0.8647 0.9179 0.9502 0.9698 0.9817 0.9889 0.9933

From the CDF function the cumulative Sigma values for r (in sigmas) equal to range from:

0-1, 0-2, 0-3, ... , 0-10

are:

0.3935, 0.6321, 0.7769, 0.8647, 0.9179, 0.9502, 0.9698, 0.9817, 0.9889, 0.9933

$\endgroup$
  • $\begingroup$ Thank you for naming the distribution. However, by (1) not differentiating between the distribution's parameter and estimates of that parameter derived from the data, (2) not stating the (strong) assumptions needed about the distribution of the balls, and (3) by being vague, you risk misleading readers. Indeed, it is unclear what the reference of your "this" is: would it be the distribution of locations of the balls? (No.) The distribution of the center of mass? (Yes, but with a scale parameter that differs from the standard deviation of the balls.) Would you like to clarify your answer? $\endgroup$ – whuber May 26 '15 at 13:45
  • $\begingroup$ filled in the gaps.... $\endgroup$ – MaxW Jun 13 '15 at 5:39
  • $\begingroup$ Thank you for the clarifications, Max. As a simple check of the correctness of your answer, let's consider one ball instead of $40$. Your answer appears to claim the distribution of the distance between this ball and the center of mass of all balls is a Rayleigh distribution. Unfortunately, in this case that distance is always zero. (The question specifically describes it as "the distance from the center of mass to each ball, which is calculated using simple geometry.") That suggests your answer may be wrong in every case, including for $40$ balls. $\endgroup$ – whuber Jun 13 '15 at 13:40
  • $\begingroup$ The distribution is about the center of mass. $\endgroup$ – MaxW Jun 13 '15 at 19:08
  • $\begingroup$ The CDF is setup for one ball of course. From the CDF 39% of the balls will fall within a circle of radius σ, 86% within 2σ, and 99% within 3σ. $\endgroup$ – MaxW Jun 15 '15 at 14:32
-1
$\begingroup$

The normal distribution, both positive and negative values, makes sense if you recognize that this normal distribution is for radius or "distance from centroid". The other variable, angle, is random and is uniformly distributed from 0-pi

$\endgroup$
  • $\begingroup$ The radius, which can never be negative, will definitely not have a Normal distribution! $\endgroup$ – whuber Sep 15 '15 at 21:10

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.