[Date Index]
[Thread Index]
[Author Index]
Re: Binomial ratio expectation
*To*: mathgroup at smc.vnet.net
*Subject*: [mg49992] Re: Binomial ratio expectation
*From*: Paul Abbott <paul at physics.uwa.edu.au>
*Date*: Mon, 9 Aug 2004 04:29:37 -0400 (EDT)
*Organization*: The University of Western Australia
*References*: <cetem8$6fs$1@smc.vnet.net>
*Sender*: owner-wri-mathgroup at wolfram.com
In article <cetem8$6fs$1 at smc.vnet.net>,
Ismo Horppu <ishorppu at NOSPAMitu.st.jyu.fi> wrote:
> I have the following problem, I need to compute
> EXPECTATION[X/(2+X)],
> where X follows Binomial distribution with n trials and success
> probability of w.
>
> I have tried to solve it with Mathematica (version 4.1) as
> Sum[((x)/(2 + x))*Binomial[n, x]*w^x*(1 - w)^(n - x), {x, 0, n}]
>
> I omit here the result which seems to be okay (according to
> simulations) for values 0<w<1. Problem is that result (intermediate
> or full simplified one) is not defined with values 0 or 1 of parameter w.
> However, it is trivial to compute the result by hand on those cases
> (as the X is then a fixed constant, 0 or n).
>
> Does anyone know how to get the full result with Mathematica, or at
> least a warning that the result is partial.
After loading the statistics stub,
<< Statistics`
the expectation value, simplified assuming 0 < w < 1 is immediate.
ev = FullSimplify[ExpectedValue[Function[x,x/(2+x)],
BinomialDistribution[n, w]], 0 < w < 1];
We need to take a limit (or use Series) as w -> 0.
Limit[ev, w -> 0]
0
However, we can just substitute in the value as w -> 1.
Simplify[ev /. w -> 1]
n/(n + 2)
> I am also interested in whether someone knows what kind of summation
> formula Mathematica uses for the sum, some kind of binomial identity
> formula perhaps? (I am unable to find which one, any references would be
> appreciated).
You can look at the code Mathematica uses to compute the ExpectedValue
for the BinomialDistribution (see Statistics`DiscreteDistributions`):
BinomialDistribution/: ExpectedValue[f_Function,
BinomialDistribution[n_, p_], opts___?OptionQ] :=
Module[{x, sum},
(
If[{opts} =!= {}, Message[ExpectedValue::sum]];
sum
) /; (sum = Sum[ Evaluate[f[x] PDF[BinomialDistribution[n, p], x]],
Evaluate[{x, 0, n}] ];
FreeQ[sum, Sum])
] /; ParameterQ[BinomialDistribution[n, p]]
If you simplify the summand you will see that you are computing
Sum[x/(x + 2) (1 - w)^(n - x) w^x Binomial[n, x], {x, 0, n}]
One way to compute this sum is to note that it can be generated by
parametric differentiation and integration from the simpler but more
general binomial sum:
FullSimplify[Sum[(1 - a)^(n - x) b^x Binomial[n, x], {x, 0, n}],
0 < w < 1]
which evaluates to
(b - a + 1)^n
Since
Assuming[x > 0, w D[1/w^2 Integrate[b b^x, {b, 0, w}], w]]
evaluates to
w^x x/(x + 2)
the sum you need to compute is just
FullSimplify[w D[1/w^2 Integrate[b (b - a + 1)^n, {b, 0, w}],w] /.
a -> w, 0 < w < 1]
and this result agrees with that from using ExpectedValue.
Cheers,
Paul
--
Paul Abbott Phone: +61 8 9380 2734
School of Physics, M013 Fax: +61 8 9380 1014
The University of Western Australia (CRICOS Provider No 00126G)
35 Stirling Highway
Crawley WA 6009 mailto:paul at physics.uwa.edu.au
AUSTRALIA http://physics.uwa.edu.au/~paul
Prev by Date:
**Re: populate a list with random numbers from normal distribution?**
Next by Date:
**Re: Re: populate a list with random numbers from normal distribution?**
Previous by thread:
**Re: Binomial ratio expectation**
Next by thread:
**Simplify or Cancel question**
| |