Re: Serious problem with CDF[MultinormalDistribution...
- To: mathgroup at smc.vnet.net
- Subject: [mg113925] Re: Serious problem with CDF[MultinormalDistribution...
- From: Myeong Ae Kang <makang at purdue.edu>
- Date: Thu, 18 Nov 2010 07:05:27 -0500 (EST)
With version 6.0 and 7.0.1, the difference is quite large. With other input
values, the results are much worse even than this.
In[1326]:= Needs["MultivariateStatistics`"]
(1.0-CDF[MultinormalDistribution[{0,0,0},Rationalize[{{1,0.559887,0.289988},
{0.559887,1,0.880157},{0.289988,0.880157,1}}]],SetPrecision
[{4.5,4.5,6.9958},50]])
(1.0-CDF[MultinormalDistribution[{0,0},{{1,0.559887},{0.559887,1}}],{4.5,4.5}])
Out[1327]= 5.00832*10^-6
Out[1328]= 6.75155*10^-6
And, similarly to this, with version 6.0, I did the calculation of CDF with 6
RV's in different orders of RV's. I just changed the orders of the sigma matrix
and range vector with zero mean. But, the results(which are suppossed to equal)
are quite different like this following:
In[1432]:= ClearAll;
Print["R,A,W ======================="]
TT={{1,0.576363,0.300013,0.888477,0.106961,0.215919},
{0.576363,1,0.888477,0.300013,0.215919,0.106961},{0.300013,0.888477,1,-
9.72249*10^-6,0.466809,-0.191847},{0.888477,0.300013,-9.72249*10^-6,1,-
0.191847,0.466809},{0.106961,0.215919,0.466809,-0.191847,1,-0.862483},
{0.215919,0.106961,-0.191847,0.466809,-0.862483,1}};
TTbnd={7.14659,7.14659,7.83599,7.83599,5.55663,5.55663};
TTp={{0,0,0,0,0,1},{0,0,0,0,1,0},{0,0,0,1,0,0},{0,0,1,0,0,0},{0,1,0,0,0,0},
{1,0,0,0,0,0}}.TT.{{0,0,0,0,0,1},{0,0,0,0,1,0},{0,0,0,1,0,0},{0,0,1,0,0,0},
{0,1,0,0,0,0},{1,0,0,0,0,0}}
Tbnd=TTbnd.Transpose[{{0,0,0,0,0,1},{0,0,0,0,1,0},{0,0,0,1,0,0},{0,0,1,0,0,0},
{0,1,0,0,0,0},{1,0,0,0,0,0}}]
1.0-CDF[MultinormalDistribution[{0,0,0,0,0,0},TTp],Tbnd]
Print["A,R,W ======================="]
TTp2={{0,0,1,0,0,0},{0,0,0,1,0,0},{0,0,0,0,1,0},{0,0,0,0,0,1},{1,0,0,0,0,0},
{0,1,0,0,0,0}}.TT.Transpose[{{0,0,1,0,0,0},{0,0,0,1,0,0},{0,0,0,0,1,0},
{0,0,0,0,0,1},{1,0,0,0,0,0},{0,1,0,0,0,0}}]
Tbnd2=TTbnd.Transpose[{{0,0,1,0,0,0},{0,0,0,1,0,0},{0,0,0,0,1,0},{0,0,0,0,0,1},
{1,0,0,0,0,0},{0,1,0,0,0,0}}]
1.0-CDF[MultinormalDistribution[{0,0,0,0,0,0},TTp2],Tbnd2]
Print["A,W,R ======================="]
TTp3={{0,0,1,0,0,0},{0,0,0,1,0,0},{1,0,0,0,0,0},{0,1,0,0,0,0},{0,0,0,0,1,0},
{0,0,0,0,0,1}}.TT.Transpose[{{0,0,1,0,0,0},{0,0,0,1,0,0},{1,0,0,0,0,0},
{0,1,0,0,0,0},{0,0,0,0,1,0},{0,0,0,0,0,1}}]
Tbnd3=TTbnd.Transpose[{{0,0,1,0,0,0},{0,0,0,1,0,0},{1,0,0,0,0,0},{0,1,0,0,0,0},
{0,0,0,0,1,0},{0,0,0,0,0,1}}]
1.0-CDF[MultinormalDistribution[{0,0,0,0,0,0},TTp3],Tbnd3]
During evaluation of In[1432]:= R,A,W =======================
Out[1436]= {{1.,-0.862483,0.466809,-0.191847,0.106961,0.215919},{-0.862483,1.,-
0.191847,0.466809,0.215919,0.106961},{0.466809,-0.191847,1.,-9.72249*10^-
6,0.300013,0.888477},{-0.191847,0.466809,-9.72249*10^-6,1.,0.888477,0.300013},
{0.106961,0.215919,0.300013,0.888477,1.,0.576363},
{0.215919,0.106961,0.888477,0.300013,0.576363,1.}}
Out[1437]= {5.55663,5.55663,7.83599,7.83599,7.14659,7.14659}
Out[1438]= 1.37516*10^-8
During evaluation of In[1432]:= A,R,W =======================
Out[1440]= {{1.,-9.72249*10^-6,0.466809,-0.191847,0.300013,0.888477},{-
9.72249*10^-6,1.,-0.191847,0.466809,0.888477,0.300013},{0.466809,-0.191847,1.,-
0.862483,0.106961,0.215919},{-0.191847,0.466809,-0.862483,1.,0.215919,0.106961},
{0.300013,0.888477,0.106961,0.215919,1.,0.576363},
{0.888477,0.300013,0.215919,0.106961,0.576363,1.}}
Out[1441]= {7.83599,7.83599,5.55663,5.55663,7.14659,7.14659}
Out[1442]= 6.50208*10^-9
During evaluation of In[1432]:= A,W,R =======================
Out[1444]= {{1.,-9.72249*10^-6,0.300013,0.888477,0.466809,-0.191847},{-
9.72249*10^-6,1.,0.888477,0.300013,-0.191847,0.466809},
{0.300013,0.888477,1.,0.576363,0.106961,0.215919},
{0.888477,0.300013,0.576363,1.,0.215919,0.106961},{0.466809,-
0.191847,0.106961,0.215919,1.,-0.862483},{-0.191847,0.466809,0.215919,0.106961,-
0.862483,1.}}
Out[1445]= {7.83599,7.83599,7.14659,7.14659,5.55663,5.55663}
Out[1446]= 1.17684*10^-14
---------------
Now, I have two questions.
1) Does mathematica version 8.0 give three similar results of the second part?
I don't have version 8.0. Before upgrading it (if possible), I wonder that.
2) Is there any other way to get the proper results with version 6.0 or 7.01?
Thank you.
Regards,
Myeong Ae Kang
=================================
Darren Glosemeyer <darreng at wolfram.com> wrote:
> On 11/17/2010 4:29 AM, Myeong Ae Kang wrote:
> > Hello.
> >
> > I have a problem, serious and urgent.
> >
> > I got wrong results of CDF of MultinormalDistribution.
> >
> > I'm calculating (1.0 -CDF[MultinormalDistribution[{0, 0,0}, {{1, 0.559887,
> > 0.289988}, {0.559887, 1,0.880157}, {0.289988, 0.880157, 1}}],
> > {4.5,4.5,6.9958}]).
> > The results is smaller than (1.0 -CDF[MultinormalDistribution[{0,0}, {{1,
> > 0.559887}, {0.559887, 1}}], {4.5,4.5}]). It should be larger.
> >
> > Actually, I have to calculate more than 6 to 12 random variables (12 x 12
> sigma
> > matrix.).
> >
> > Mathematica 6.0 and 7.0.1 gave me the same wrong result.
> >
> > I compared this with the result of another system.
> > That system gave the expected value.
> >
> > To fix this problem, what should I do?
> >
> >
> > Regards,
> > Myeong Ae Kang
> >
> >
> >
>
>
> This is addressed in version 8. Here are the results:
>
>
>
> In[1]:= threeD = (1.0 -
> CDF[MultinormalDistribution[{0, 0,
> 0}, {{1, 0.559887, 0.289988}, {0.559887, 1,
> 0.880157}, {0.289988, 0.880157, 1}}], {4.5, 4.5, 6.9958}])
>
> -6
> Out[1]= 6.75155 10
>
> In[2]:= twoD = (1.0 -
> CDF[MultinormalDistribution[{0,
> 0}, {{1, 0.559887}, {0.559887, 1}}], {4.5, 4.5}])
>
> -6
> Out[2]= 6.75155 10
>
> In[3]:= threeD > twoD
>
> Out[3]= True
>
>
> The difference is at the level of machine precision numerical error:
>
>
> In[4]:= threeD - twoD
>
> -16
> Out[4]= 1.11022 10
>
>
>
> Note that if we artificially increase the precision as in the example
> below, we get a difference at the same order of magnitude, so the
> precision in the difference above is as good as could be expected (and
> in fact probably better than should be expected) at machine precision.
>
>
> In[5]:= CDF[MultinormalDistribution[{0, 0, 0},
> Rationalize[{{1, 0.559887, 0.289988}, {0.559887, 1,
> 0.880157}, {0.289988, 0.880157, 1}}, 0]],
> SetPrecision[{4.5, 4.5, 6.9958}, 50]]
>
> Out[5]= 0.999993248446123977600817893853727681422304494350
>
> In[6]:= CDF[MultinormalDistribution[{0, 0},
> Rationalize[{{1, 0.559887}, {0.559887, 1}}, 0]],
> SetPrecision[{4.5, 4.5}, 50]]
>
> Out[6]= 0.999993248446124136306526964076086153533229885505
>
> In[7]:= % - %%
>
> -16
> Out[7]= 1.5870570907022235847211092539116 10
>
>
> Darren Glosemeyer
> Wolfram Research
>