Uncertainty for a star being between two values

[omoplata]

I have a set a of stars, each of which have a magnitude in the K band ( K ) and a 1 sigma uncertainty ( dK ) associated with it. For a program I'm writing, I need to see whether each star lies within a maximum K magnitude (K_max ) and a minimum K magnitude ( K_min ). I also need to find a probability ( p ) for this star to lie between these two values.

For example, if star 1 has magnitude K_1, I need to find what is the probability p for K_1 to lie between K_min and K_max. I have the 1 sigma uncertainty dK_1 for this magnitude.


Extra info:

The reason I want this probability p is to use it as an uncertainty. For each star, I'm going to check whether it's K lies within K_max and K_min, and assign a value F = 1 if it does, and F = 0 if it does not. Then I want to use p to find the uncertainty of F, or dF. If F = 1, then I'll use dF = p. If F = 0, then I'll use dF = (1-p). Why I want dF is to use it in an adaptive binning program, that someone has already written.

 

[Stephen Tashi]

One fundamental question is (for the sake of us non-astronomers) is whether the K band magnitude is a constant in time or at least in whatever time period these programs are concerned with. In other words, what physical units does the standard deviation data have? Is it in units of energy? Or is a unit of time involved? I don't understand whether the fact that the magnitude has a standard deviation indicates an uncertainty in measurement or whether it indicates an actual physical variation. Can you give us a good link to a explanation of K band magnitudes? Is it known that these mangitudes are normally distributed?

 

[omoplata]

One fundamental question is (for the sake of us non-astronomers) is whether the K band magnitude is a constant in time or at least in whatever time period these programs are concerned with. In other words, what physical units does the standard deviation data have? Is it in units of energy? Or is a unit of time involved? I don't understand whether the fact that the magnitude has a standard deviation indicates an uncertainty in measurement or whether it indicates an actual physical variation. Can you give us a good link to a explanation of K band magnitudes? Is it known that these mangitudes are normally distributed?

Thank you for the interest in my problem.

The K band magnitude I speak about is the apparent magnitude, only within a certain wavelength range. The apparent magnitude is defined to be,

$$m_x = -2.5 \cdot \log_{10} \left( \frac{F_x}{F_x^0} \right)$$
where $F_x$ is the flux of the star in the band x (flux received by the telescope), and $F_x^0$ is some reference flux also in the band x, the value of which is known, most commonly that of Vega. In this case, the band is K. But in the initial statement of my problem, I named $m_k$ as K.

Here's the wikipedia link for apparent magnitude( http://en.wikipedia.org/wiki/Apparent_magnitude )

K band is just when we apply a filter to the telescope so it only receives the light between approximately 2000 nm and 2400 nm. The relevant wikipedia link is ( http://en.wikipedia.org/wiki/Photometric_system ).

But I don't think the above is really important for this problem. I think for this problem we can consider K to be a random variable which can take any real value.

So, it doesn't have units. And for the time period the program is concerned with, it is constant w.r.t. time.

The uncertainty I gave is the uncertainty of measurement.

The brighter a star is seen by the observer, the lower its numerical value of the apparent magnitude is ( goes on to negative infinity ). The dimmer a star is seen by an observer, the higher its numerical value of the apparent magnitude is ( goes on to positive infinity ). As the apparent magnitude becomes more negative, there are less and less stars ( The sun is the brightest star we see ). As the apparent magnitude becomes more positive, the number count of stars increase. There are many more dimmer stars. The highest apparent magnitude ( dimmest stars ) we can see are limited by the sensitivity of our telescopes. So I guess that means the K apparent magnitudes are NOT normally distributed.

Assuming the stars are uniformly distributed throughout space ( they are clustered together in galaxies, but assuming we are counting stars within a small volume in our galaxy ) the CUMULATIVE number count of stars can be given to be $N = c \cdot 10^{3 K/5}$, where c is some constant. The derivation of that is given in this tread < https://www.physicsforums.com/showthread.php?t=480112 > ( Please note that the notation is different. Also it looks weird because PF LaTeX was changed since that tread was made). N is the number of stars counted up to apparent magnitude K from apparent magnitude negative infinity. ( $N = \int_{L = -\infty}^{K} dN$ where dN is the number of stars between apparent magnitude L and L + dL )

 

[Stephen Tashi]

From those links and reading a few pages available on Google books from "Observing The Universe" edited by Andrew Norton, I conclude (from Norton) that the uncertainty in measurement of magnitudes of a particular star is treated as a normal distribution (even though the distribution of magnitudes in the population isn't). Norton doesn't give a mathematical treatment but his rules of thumb (e.g. that plus or minus one uncertainty corresponds to roughly 2/3 of the observed values) suggest that his book is doing that.

If you want to (literally) compute the probability that the actual magnitude of a particular star is between a given $K_{min}$ and $K_{max}$ then you must assume a scenario that acknowledges that the value of the actual magnitude is probabilistic. This would be a Bayesian statistical approach. You can't use the typical scenario from "frequentist" statistics where the actual magnitude is a "fixed but unknown value".

I'm a fan of Bayesian statistics, however it would surprise me if the binning program really wanted a probability based on a Bayesian statistical approach. Given the common use of "frequentist" statistics, it is more likely that it wants the probability under the standard normal density curve ( mean 0 ,sigma = 1) associated with "confidence interval" of $[\frac{K_{min} - K_{obs}}{\sigma}, \frac{K_{max}-K_{obs}}{\sigma}]$.

(The terminology for a superficially similar interval in Bayesian statistics is a "credible interval". )

In browsing the web, I was surprised to find that astronomy uses very sophisticated statistical methods such as the "Malmquist correction", so I may be underestimating the binning program. Nevertheless, I wouldn't trust how the program documentation defined "p" unless it throroughly explained how "p" is used.

To do the Bayesian calculation, you would define a problem like this:

Given this particular star was picked at random from the population of stars and given the observed magnitude was $K_{obs}$ and that the measurement errors are normally distributed with standard deviation $\sigma$, what is the probability that the actual magnitude is in the interval $[ K_{min}, K_{max}]$. This calculation can be done numerically. We can attempt to do the details if this is really the number you want.

 

[omoplata]

Given this particular star was picked at random from the population of stars and given the observed magnitude was $K_{obs}$ and that the measurement errors are normally distributed with standard deviation $\sigma$, what is the probability that the actual magnitude is in the interval $[ K_{min}, K_{max}]$. This calculation can be done numerically. We can attempt to do the details if this is really the number you want.

Yeah, I think that's what I want.

I guess I can find the probability by assuming a normal distribution with mean $K_{obs}$ and standard deviation $\sigma = dK_{obs}$, and then integrate from $K_{min}$ to $K_{max}$?

 

[Stephen Tashi]

Yeah, I think that's what I want.

I guess I can find the probability by assuming a normal distribution with mean $K_{obs}$ and standard deviation $\sigma = dK_{obs}$, and then integrate from $K_{min}$ to $K_{max}$?

Doing that doesn't compute what you said you wanted. I think that numerically, it might give you approximately what you wanted. (I also think that this may be what the binning program expects, but you must verify that.)

If you do that integration and interpret it as a statement about the probability of the actual magnitude being in the interval, you'll have trouble justifying why you did it. Do you know that the mean of the distribution of the actual magnitude for that star is equal to $K_{obs}$? That's what you'd be assuming.

There are two ways that the integration you suggested can pop up in statistics. It pops up in "frequentist" statistics as a calculation of of "confidence" but the result can't be interpreted as a statement about the probability of the actual magnitude being in a specific interval, an interval that has numerical endpoints.

It also can occur in a Bayesian approach where you assume something resembling complete ignorance about the magnitudes of stars. Of course, you can't make a probability distribution that says "the magnitude is equally likely to be any real number". You can pick some large number L and assume the magnitude is equally likely to be between -L and L. Then you can calculate the probability distribution for the magnitude given that the magnitude you observed was $K_{obs}$. Then integrate that distribution from $K_{min}$ to $K_{max}$. It comes out to be approximately the same as calculation that you suggested.

The full Bayesian calculation may be too complicated for you. It has the following pattern:

Let F(x) be the "prior" cumulative distribution of star magnitudes. Let f(x) be its density.

For F(x) you could use the distribution in other thread you mentioned. Or you could assume the distribution of magnitudes is uniform on -L to L. Or you could even assume that "the population" will consist of some stars that you judge to be similar to the particular star in question and you could use the distribution of their magnitudes as F(x).

Let n(x,a,sigma(a)) denote a normal density function with mean 'a' and standard deviation,sigma. I have made sigma a function of 'a'. If things are simple then the standard deviation of the measurement error would be the same for all magnitudes of measurements. I suspect it isn't. You're the expert in that question.

To compute the conditional density g(y) for the actual magnitude y given the observation $K_{obs}$, I think of it as

g(y) = Pr( actual magnitude = y | observed magnitude = $K_{obs}$)

= Pr( actual magnitude = y and observed magnitude = $K_{obs}$) / Pr( observed magnitude = $K_{obs}$

= n( $K_{obs}$, y, sigma(y))/ $\int$ f(x) n($K_{obs}$, x, sigma(x)) dx

To find the probability that the actual magnitude is in [$K_{min} , K_{max}$ ] compute $\int_{K_{min}}^{K_{max}} g(y) dy$.

 

[omoplata]

Hi. Sorry about the late reply. I had some urgent stuff to take care of. I really appreciate your help this problem.

Let F(x) be the "prior" cumulative distribution of star magnitudes. Let f(x) be its density.

For F(x) you could use the distribution in other thread you mentioned. Or you could assume the distribution of magnitudes is uniform on -L to L. Or you could even assume that "the population" will consist of some stars that you judge to be similar to the particular star in question and you could use the distribution of their magnitudes as F(x).

There's a problem with the distribution that I mentioned in the other thread. It's maximum value is not 1. It goes on to infinity. So it cannot be used as a cumulative distribution function. So does its derivative, so that can't be used as the density function. Is there any way of modifying that distribution so that it's suitable for use?

To compute the conditional density g(y) for the actual magnitude y given the observation $K_{obs}$, I think of it as

g(y) = Pr( actual magnitude = y | observed magnitude = $K_{obs}$)

= Pr( actual magnitude = y and observed magnitude = $K_{obs}$) / Pr( observed magnitude = $K_{obs}$

= n( $K_{obs}$, y, sigma(y))/ $\int$ f(x) n($K_{obs}$, x, sigma(x)) dx

To find the probability that the actual magnitude is in [$K_{min} , K_{max}$ ] compute $\int_{K_{min}}^{K_{max}} g(y) dy$.

I get the first two lines. But I don't understand the following.

Why does Pr( actual magnitude = y and observed magnitude = $K_{obs}$) = n( $K_{obs}$, y, sigma(y)) ?

Why does Pr( observed magnitude = $K_{obs}$) = $\int$ f(x) n($K_{obs}$, x, sigma(x)) dx ?

 

[Stephen Tashi]

There's a problem with the distribution that I mentioned in the other thread. It's maximum value is not 1. It goes on to infinity. So it cannot be used as a cumulative distribution function. So does its derivative, so that can't be used as the density function. Is there any way of modifying that distribution so that it's suitable for use?

I see. It's a "count", not a distribution. If you can bound the population of stars you are dealing with some maximum magnitude, you could divide the count by the count of stars below that magnitude to get a probability distribution. If you are willing to assume there is some maximum possible error in an measurement of magnitude $K_{obs}$, you can consider only the population stars within those bounds. You can have a different F(x) for each $K_{obs}$.

I get the first two lines. But I don't understand the following.

Why does Pr( actual magnitude = y and observed magnitude = $K_{obs}$) = n( $K_{obs}$, y, sigma(y)) ?

You're correct not to understand it because I made a mistake!

It should be f(y) n($K_{obs}$, y, sigma(y)). It's just like the denominator except we aren't integrating over the variable y.

Then the expression represents Pr( true magnitude = y) Pr( we observe magnitude $K_{obs}$ given that the true magnitude is y).
So its an example of P(A and B) = P(A) P(B given A)

Why does Pr( observed magnitude = $K_{obs}$) = $\int$ f(x) n($K_{obs}$, x, sigma(x)) dx ?

One way of interpreting the Bayesian view is that there are many possible worlds. In each different world, there is a different true magnitude x of the star. When the magnitude is measured, the error in measurement produces a value that is normally distributed with mean x and variance sigma(x). The "probability of observing $K_obs$" refers to the probability of observing it in the following process: Pick a world at random according to the density f(x). Then make an observation of the stars magnitude in that world. The probability of observing $K_obs$ is, in a manner of speaking, the sum of the probabilities of observing it taken over all worlds. So it is the sum of (probability of being in the world where the x and sigma(x) define the true distribution of measurement errors) ( probability of observing $K_{obs}$ given that we are in that world). The integrand has the same pattern P(A) P(B| A) as the corrected numerator.

 

[omoplata]

It should be f(y) n($K_{obs}$, y, sigma(y)). It's just like the denominator except we aren't integrating over the variable y.

Then the expression represents Pr( true magnitude = y) Pr( we observe magnitude $K_{obs}$ given that the true magnitude is y).
So its an example of P(A and B) = P(A) P(B given A)

I have a question about calculating probability here. My understanding is the the probability distribution function cannot be directly used to calculate the probability for a random variable to be at a certain value. To find the probability for a random variable to be between two values, the p.d.f. has to be integrated between the two values. Am i wrong?

For example, if h(x) is a p.d.f., to find the probability to be between x_1 and x_2 we have to integrate h(x) from x_1 to x_2 w.r.t. x.

But here, P(A) is directly taken to be f(y), and P(B|A) is directly taken to be n($K_{obs}$, y, sigma(y)), without integrating between any y_1 and y_2 values.

One way of interpreting the Bayesian view is that there are many possible worlds. In each different world, there is a different true magnitude x of the star. When the magnitude is measured, the error in measurement produces a value that is normally distributed with mean x and variance sigma(x). The "probability of observing $K_obs$" refers to the probability of observing it in the following process: Pick a world at random according to the density f(x). Then make an observation of the stars magnitude in that world. The probability of observing $K_obs$ is, in a manner of speaking, the sum of the probabilities of observing it taken over all worlds. So it is the sum of (probability of being in the world where the x and sigma(x) define the true distribution of measurement errors) ( probability of observing $K_{obs}$ given that we are in that world). The integrand has the same pattern P(A) P(B| A) as the corrected numerator.

OK. So in,

Pr( observed magnitude = Kobs) = ∫ f(x) n(K_obs, x, sigma(x)) dx ,

f(x) corresponds to P(A), and n(K_obs, x, sigma(x)) corresponds to P(B|A)? P(B|A) is a normal distribution with mean x and standard deviation sigma(x) because in that particular world, the mean is x and and x values are normally distributed with standard deviation sigma(x)? So A means the actual value is x, while B means the observed value is K_obs?

 

[Stephen Tashi]

I have a question about calculating probability here. My understanding is the the probability distribution function cannot be directly used to calculate the probability for a random variable to be at a certain value. To find the probability for a random variable to be between two values, the p.d.f. has to be integrated between the two values. Am i wrong?

You are correct that saying the value of a continuous density f(x) as "The probability of x" is not correct. For example, the value of the density of the uniform distribution on [0, 1/2] must be f(x) = 2. Nevertheless, mis-thinking of f(x) as "the probability of x" usually leads to the correct expressions in most probability calculations. It's somewhat like "illegal" manipulations of the symbols dy and dx in calculus. They help you remember the correct results.

For example, if h(x) is a p.d.f., to find the probability to be between x_1 and x_2 we have to integrate h(x) from x_1 to x_2 w.r.t. x.


But here, P(A) is directly taken to be f(y), and P(B|A) is directly taken to be n($K_{obs}$, y, sigma(y)), without integrating between any y_1 and y_2 values.

Your objection is correct. Getting the result in the rigourous manner would be an interesting challenge. I think there are results from the "measure theory" that would justify the calculation, but although I took that course years ago, I can't quote them. At a lesser level of rigor, we could look at computing things like Pr( we observe a magnitude in [K_{obs{ -dk, K_{obs} + dk}] and the true magnitude is in [y - dy, y + dy] and try to reason out the answers in those terms.

OK. So in,

Pr( observed magnitude = Kobs) = ∫ f(x) n(K_obs, x, sigma(x)) dx ,

f(x) corresponds to P(A), and n(K_obs, x, sigma(x)) corresponds to P(B|A)? P(B|A) is a normal distribution with mean x and standard deviation sigma(x) because in that particular world, the mean is x and and x values are normally distributed with standard deviation sigma(x)? So A means the actual value is x, while B means the observed value is K_obs?

Yes

------

Incidentally, there is an interesting book by the late E.T. Janes (who I would characterize as a flaming Bayesian) called "Probablity Theory: The Logic Of Science". It's available as a PDF at various places on the web and its is written from the viewpoint of a physicist.