Re: Strange statistics function integration
- To: mathgroup at smc.vnet.net
- Subject: [mg59600] Re: Strange statistics function integration
- From: "antononcube" <antononcube at gmail.com>
- Date: Sat, 13 Aug 2005 03:26:35 -0400 (EDT)
- References: <ddh7v0$b7e$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
If PosterioriWithNormalProb1 is called with exact numbers we get a
correct result:
In[32]:= PosterioriWithNormalProb1[1/2,1/1000][8,0]//N
Out[32]= 0.00390669
Since PosterioriWithNormalProb1 is called with numerical values it is a
good idea, at least in general, to use NIntegrate instead of Integrate.
NIntegrate though indicates loss of precision:
In[1]:= << Statistics`ContinuousDistributions`
In[2]:=PosterioriWithNormalProb1N[mu_?NumberQ,sigma_?NumberQ][n_?NumberQ,k_?NumberQ]:=
Binomial[n, k]*
NIntegrate[
x^k*(1 - x)^(n - k)*
PDF[NormalDistribution[mu, sigma], x], {x, -Infinity,
Infinity}];
In[3]:=PosterioriWithNormalProb1N[0.5, 0.001][8, 0]
NIntegrate::ploss:
Numerical integration stopping due to loss of precision. Achieved
neither the requested PrecisionGoal nor AccuracyGoal; suspect one of
the following: highly oscillatory integrand or the true value of the
integral is 0. If your integrand is oscillatory on a (semi-)infinite
interval try using the option Method->Oscillatory in NIntegrate.
Out[3]= 0.
This is due to the fact that the integrand is basically a very narrow
spike, and it gets undersampled. Increasing MinRecursion and
MaxRecursion gives a result that agrees with the Integrate result with
exact numbers above:
In[4]:= PosterioriWithNormalProb1N[mu_?NumberQ,
sigma_?NumberQ][n_?NumberQ, k_?NumberQ] :=Binomial[n,
k]*NIntegrate[x^k*(1 - x)^(n - k)*PDF[NormalDistribution[mu, sigma],
x], {x, -Infinity, Infinity},MinRecursion -> 2, MaxRecursion -> 12];
In[5]:=PosterioriWithNormalProb1N[0.5, 0.001][8, 0]
NIntegrate::slwcon:
Numerical integration converging too slowly; suspect one of the
following: singularity, value of the integration being 0, oscillatory
integrand, or insufficient WorkingPrecision. If your integrand is
oscillatory try using the option Method->Oscillatory in NIntegrate.
Out[5]= 0.00390669