Power calculations for 2 by m tables in one and two stage

Download Report

Transcript Power calculations for 2 by m tables in one and two stage

Bounds on the probability
of the union of events
International Colloquium on Stochastic Modeling and Optimization
dedicated to the 80th Birthday of Professor András Prékopa
József Bukszár
Medical College of Virginia, Virginia Commonwealth University
email: [email protected]
The underlying problem
Let  A1 ,..., An  be arbitrary events.
Our purpose is to give lower and upper bounds for
Pr A1 ... An 


based on intersection probabilities ...,Pr Ai1 ... Aik ,...,
where
k
an integer given in advance.
Example: second and first order Bonferroni-bounds
 Pr A    PrA  A   Pr A  ... A    Pr A 
n
i 1
n
i
1i  j  n
i
j
S2
1
n
i 1
i
S1
Motivation
Estimate values of multivariate distribution functions:
F x1,...,xn   Pr X1  x1 ,..., X n  xn   Pr A1  ... An  
A1

An
1  Pr A1c  ... Anc

We can give lower and upper bounds for F x1 ,...,xn 
based on the marginal distribution function values.
For that we need tight bounds that are based on only few
intersection probabilities.
Motivation
Network reliability
Vertices represent stations.
Edges represent phone lines,
each of them is busy with a
certain probability.
X
U
Y
What is the probability that
we can call Y from X?
V
A1 is the event that we can call Y from X through the red line.
The probability that we can call Y from X is Pr(A1 U…U An)
Each event corresponds to a path connecting X and Y.
Hunter-Worsley boundA
Given a complete graph whose
vertices are identified by the events.
A2
1
A3
A7
The weight of the edge connecting
Ai and Aj is Pr(Ai  Aj ).
Let T = (V,E) be the maximum weight
spanning tree on the complete graph.
Hunter-Worsley bound:
n
Pr A1  ... An    Pr Ai  
i 1
A4
A6
A5
Pr(A4  A5 )
PrA  A 

 
i , j E
i
j
edges of the maximum
weight spanning tree
Advantages:


quickly computable: the maximum weight spanning tree
can be obtained by a greedy algorithm (e.g. by Prim’s
algorithm)
requires only n-1 intersection probabilities
Disadvantages:


no improvement is available (problematic if the bound is
greater than 1)
provides upper bound only
m-multicherry
v1
m > 0 integer
DEF: An m-multicherry is a
hypergraph of the form
(V,E2,…,Em+1),
where
m - multicherry
v1 ,...,vm , vm1 
v2
vm+1 (middle
vertex)
vm-1
V=(v1,…,vm+1) is the set of vertices,
vm
Ei = { H | vm+1 H  {v1,…,vm}, |H|=i }
are the set of hyperedges.
Example: 2-multicherry (cherry)
v1
E2 = { {v1,v3}, {v2,v3} }
E3 = { {v1,v2,v3} }
v2
v3
m-multitree (recursive def.)
DEF: An m-multitree is a hypergraph of the form
(V,E2,…, Ei,…,Em+1),
set of vertices
set of hyperedges with i vertices
i.)The smallest m-multitree has m vertices.
Ei is the set of all subsets of V with i vertices (Em+1 = ).
ii.) From an m-multitree we can obtain another m-multitree by
adding a new vertex v and an m-multicherry with middle vertex v.
”old multitree”
v
Example of a 2-multitree (cherry tree)
Building up a cherry tree (V,E2,E3):
1
1
2
1
3
2
2
E2 = { {1,2}, {1,3}, {2,3} }
E3 = { {1,2,3} }
E2 = { {1,2} }
E3 = 
2
5
6
3
3
1
4
4
E2 = { {1,2}, {1,3}, {2,3}, {2,4}, {3,4} }
E3 = { {1,2,3}, {2,3,4} }
5
1
3
2
4
7
Recursion is not unique.
We can add the vertices
for example in the order
2,3,7,4,1,5,6 to obtain
the same cherry tree.
Upper bounds by m-multitrees
DEF: The weight of an m-multitree =(V,E2, …,Em+1) with V={1,…,n}
is defined by
w 
PrA

 
i1
i1 ,i2 E2
...   1

 Ai2 
m 1
 PrA
i1 ,i2 ,i3 E3
 Pr A
i1 ,...,im1 Em1
i1
i1

 Ai2  Ai3  ...

 ...  Aim1 .
THEOREM: For any m-multitree =(V,E2, …,Em+1) with V={1,…,n}
we have
n
Pr A1  ... An    Pr Ai   w   S1  w .
i 1
Special case m=1 provides us the Hunter-Worsley bound.
Some properties of m-multitrees
m-multitrees provide us m+1 order upper bounds.
An m-multitree is completely determined by its set of vertices
and set of edges. In other words, although an m-multitree is a
hypergraph, it can be identified by its graph.
A6
A1
A2
A5
A3
A4
A7
There are only O(n) intersection probabilities involved in an
m-multitree bound. The number of intersection probabilities with
at most m+1 events is O(nm+1), an m-multitree bound uses
only O(n) out of them. Useful when the intersection probabilities
have to be evaluated, e.g. estimate of multivariate distr. function.
Unfortunately, the greedy algorithm generally does not provide
the maximum weight m-multitree if m > 1.
THEOREM: An m-multitree =(V,E2, …,Em+1) with V={1,…,n}
can be extended to an m+1-multitree ’=(V,E’2, …,E’m+2) with
w(’) w().
(Extension means that Ei  E’i for every i.)
ALGORITHM TO FIND HEAVY M-MULTITREE
The theorem enables us to find a heavy m-multitree by starting
from the heaviest (1-multi)tree and increasing m step by step
to improve on the bound.
Only the intersection probabilities involved in the multitree bound
are needed, i.e. O(n).
Lower bounds for a 30-variate normal distribution function value
bound
sec.
Sobel-Uppuluri-Galambos
-1.949941
45.59
2-mutlitree
0.831257
0.05
4-matroid tree (Grable)
0.550906
46.41
3-mutlitree
0.849664
0.11
S1,S2 sharp
0.605092
0.03
4-mutlitree
0.857520
0.45
S1,S2,S3 sharp
0.750745
2.09
5-mutlitree
0.861284
2.19
S1,S2,S3,S4 sharp
0.792341
47.95
6-mutlitree
0.863209
4.89
S1,S2,S3 multitree aggregated
0.681811
2.02
7-mutlitree
0.864234
10.60
S1,S2,S3,S4 multitree aggregat.
0.728691
45.59
8-mutlitree
0.864794
36.58
Hunter-Worsley (1-multitree)
0.776914
0.03
UPPER B.
0.865619
90.12
name
Sk 
 PA
1i1 ...ik n
i1
 ... Aik

name
bound sec.
Marginal function values were computed
by Genz’s Fortran code SADMVN and
IMSL subroutines MDNOR and MDBNOR.
x1=1.55, x2=1.6, …,x29=2.95, x30=3.0
0.8 if i  j
Covariance rij =
1.0 if i = j
Lower bounds for a 30-variate normal distribution function value
bound
sec.
Sobel-Uppuluri-Galambos
-0.062999
48.23
2-mutlitree
0.943375
0.06
4-matroid tree (Grable)
0.808681
58.93
3-mutlitree
0.947100
0.09
S1,S2 sharp
0.866069
0.03
4-mutlitree
0.949233
0.23
S1,S2,S3 sharp
0.918326
2.47
5-mutlitree
0.950598
1.54
S1,S2,S3,S4 sharp
0.930921
50.65
6-mutlitree
0.951570
4.27
S1,S2,S3 multitree aggregated
0.893123
2.37
7-mutlitree
0.952281
18.69
S1,S2,S3,S4 multitree aggregat.
0.909447
48.27
8-mutlitree
0.952823
77.38
Hunter-Worsley (1-multitree)
0.934630
0.03
UPPER B.
0.955939
31.04
name
x1 = x2 = … = x29 = x30= 2.5
Covariance rij 
name
bound sec.
Normal random variables with this
covariance matrix are used for
estimating American Option Price.
i j for all 1 i  j  30.
(h,m)-hypermultitree (recursive def.)
h  0, m  1 integers
DEF: An (h,m)-hypermultitree is a hypergraph of the form
(V, hE2,…, hEi,…,hEm+1),
set of vertices
set of hyperedges with h+i vertices
i.) (0,m)-hypermultitrees are the same as m-multitrees, 0Ei  Ei .
ii.) The smallest (h,m)-hypermultitree has h+m vertices, and
hEi is the set of all subsets of V with h+i vertices (hEm+1 = ).
iii.) From an (h,m)-hypermultitree =(V, hE2,…,hEm+1), we can obtain
another (h,m)-hypermultitree ’=(V’, hE’2,…,hE’m+1),
by adding a new vertex v and some hyperedges in the following way.
Let =(V, h-1E*2,…,h-1E*m+1) be an arbitrary (h-1,m)-hypermultitree (on ).
We add the hyperedges of  extended by v to obtain ’, i.e.
V '  V  v

'
*


E

E

H

v
|
H

E
h
i
h i
h 1
i
and

Example:
(1,1)-hypermultitree
=( {a,b,c,d}, 1E2)
”old multitree”
Take a (0,1)-hypermultitree (i.e. tree)
=( {a,b,c,d}, { {a,b},{a,c},{a,d} } ) on .
b
a
c
d
Hyperedges added at this step:
{a,b,v}, {a,c,v} and {a,d,v}.
v
Bounds by (h,m)-hypermultitrees
DEF: The weight of an (h,m)-hypermultitree =(V,hE2, …,hEm+1)
with V = {1,…,n} is defined by
w 
 PrA
i1 ,...,ih2
 h E2
i1

 ... Aih2 
...   1
m 1
 PrA
i1 ,...,ih3
 Pr A
i1 ,...,ihm1
 h E m1
i1
 h E3
i1

 ... Aih3  ...

 ...  Aihm1 .
THEOREM: For any (h,m)-hypermultitree =(V, hE2, …, hEm+1)
with V = {1,…,n} the following inequalities hold
i.) if h is even Pr A1  ... An   S1  S2  ...  Sh1  w,
ii.) if h is odd Pr A1  ... An   S1  S2  ...  Sh1  w,
where Sk 
 PA
1i1 ...ik n
i1

 ... Aik .
Special case m = 1 Tomescu bounds; h = 0 the multitree bounds.
Some properties of (h,m)-hypermultitrees




h+m+1 order bounds,
lower bounds if h is even, upper bounds if h is odd,
the heavier is the hypermultitree, the better is the
bound,
based on O(nh+1) intersection probabilities.
Remark: Consequently, for upper bounds h = 0,
for lower bounds h = 1 is a cost-effective choice,
especially in applications where the intersection
probabilities have to be evaluated.
THEOREM: An (h,m)-hypermultitree =(V, hE2 , …, hEm+1)
with V={1,…,n} can be extended to an (h,m+1)-hypermultitree
’=(V, hE’2, …,hE’m+1) with w(’) w().
(Extension means that Ei  E’i for every i.)
ALGORITHM TO FIND HEAVY (1,m)-HYPERMULTITREE
We find a heavy (1,1)-hypermultitree by a greedy algorithm.
Based on the theorem we extend this (1,1)-hypermultitree to
a (1,2)-hypermultitree that we extend to a (1,3)-hypermultitree etc.
At the end of the algorithm we obtain a (1,m)-hypermultitree.
This stepwise extension can be done in a single step, i.e.
the initial (1,1)-hypermultitree can be extended in a single step to
a (1,m)-hypermultitree that has higher weight.
Short formulae for the bounds
There is a short formula to compute m-multitree
( (1,m)-hypermultitree ) bounds containing n-m
( n-m
2 ) intersection probabilities of the type


Pr Ai1 ... Aik  A ... A .
c
j1
c
jl
altogether m+1 (m+2) events
In other words, there are some complement events
included in the above intersection probabilities.
They can be evaluated in applications where bounds
for values of multivariate distribution function values
are sought.
Lower bounds
Upper bounds
seconds
seconds
4-matroid tree (Grable)
428.86
0.719174
0.972666
13.40
3-matroid tree (Grable)
Hunter-Worsley (1-mul)
0.01
0.861747
2-multitree
0.05
0.877985
0.972666
1.14
Tomescu ( (1,1)-hyperm.)
3-multitree
0.07
0.884730
0.909455
1.27
(1,2)-hypermultitree
4-multitree
0.08
0.888459
0.903721
1.43
(1,3)-hypermultitree
5-multitree
0.20
0.890727
0.900948
1.65
(1,4)-hypermultitree
6-multitree
0.77
0.892263
0.899480
2.46
(1,5)-hypermultitree
7-multitree
2.04
0.893355
0.898524
3.08
(1,6)-hypermultitree
8-multitree
8.63
0.894148
0.897911
5.43
(1,7)-hypermultitree
x1=1.84, x2=1.88, …,x29=2.96, x30=3.0
Covariance
rij 
for all 1 i  j  30.
i j
Computation were made by a
CELERON II 850MHz computer.
Marginal function values were computed
by Genz’s Fortran code SADMVN and
IMSL subroutines MDNOR and MDBNOR.
Simulating multivariate normal
distribution function values
Tamás Szántai developed and implemented a method to simulate
multivariate normal distribution function values based on multitrees
and hypermultitrees.
The code simulates the difference between a lower (upper) bound
and the real function value and calculates  () / see the figure /.
Szántai showed that  and  are negatively correlated unbiased
estimators, thus
a + b
is an unbiased estimator of the function value with lower variance,
where a+b =1, a>0,b>0. Values of a and b are chosen optimally
(variance is minimized).
real value
 = simulation based
on lower bounds
 = simulation based
on upper bounds
Simulating multivariate normal
distribution function values (cont’d)
Another version: Let  be simulated function value obtained by
the crude Monte-Carlo method. Then
a + b + c 
is an unbiased estimator of the function value, where a+b+c =1,
a > 0, b > 0 and c > 0.
Szántai’s code turned out to be several thousand times more
effective than the crude Monte-Carlo simulation when the function
value is high and the dimension is 20-50.
The gain in effectiveness is somewhat less but still significant
for medium (low) function values 20-50 (20-30) dimension.
Some care must be taken to select m for the m-multitree
( (1,m)-hypermultitree ) bounds.
t-cherry trees
4
2
3
5
1
t-cherry tree
not t-cherry tree
vertex 1 and 4 are not adjacent
DEF: A cherry tree (2-multitree) is called a t-cherry tree if the
two non-middle vertices of every cherry are adjacent.
t-cherry trees (cont’d)
THEOREM: A t-cherry tree bound can always be identified as
the objective function value of the dual feasible basis in the
Boolean probability bounding problem.
REM: The Boolean probability bounding problem is a linear
programming problem with 2n - 1 number of variables
(n is the number of events).
REM: The same is not true for an arbitrary cherry tree.
CONJECTURE: The above theorem can be generalized to
m-multitrees.
Open Questions
Are there tight lower bounds for Pr(A1… An) of arbitrary order
that are based on O(n) number of intersection probabilities?
Is there a polynomial time algorithm that finds the maximum weight
m-multitree if m > 1? If not, then can the family of all m-multitrees
on n vertices be extended to a matroid?
Same question for (h,m)-hypermultitrees.
What are the best lower or upper bounds of a certain order?
We have seen that t-cherry trees provide us the best third order
upper bounds on certain examples, but not on all of them.
The underlying Stoch. Optim. Problem
h( x)  min
Subject to
h0 ( x)  Prg1 ( x, Y )  0,...,gn ( x, Y )  0  p
h1 ( x)  p1 ,...,hm ( x)  pm ,
Where Y is a random variable with known distribution, and
p is a constant, typically between 0.9 and 1.
This is the probability of the intersection of events gi ( x, Y )  0 .
Applying lower (upper) bound instead of the intersection probability
shrinks (extends) the set of feasible solutions.
Strategy 1 (based on lower bounds):
1.
2.
Solve the problem with a lower bound in the place of the
intersection probability
Iterate Step 1. using a better bound until optimality holds or
using the original probabilistic constraint
Strategy 2 (based on upper bounds):
1.
2.
Solve the problem with an upper bound in the place of the
intersection probability
Iterate Step 1. using a better bound until feasibility holds or
using the original probabilistic constraint
As another application, Tamás Szántai restricted the search interval with bounds
in his line search method to find the boundary points of feasible solutions.
Prékopa’s theorem
If g1 ( x, y),...,g n ( x, y) are concave functions and Y has a
continuous probability distribution with logarithmically concave
probability density function, then the function
h0 ( x)  Prg1 ( x, Y )  0,..., gn ( x, Y )  0
is also logarithmically concave.
Corr.: the set of x satisfying the probabilistic constraints
h0 ( x)  Prg1 ( x, Y )  0,...,gn ( x, Y )  0  p
is convex.
Exa: Let the probabilistic constraints in the underlying problem be
PrTi x  Yi , i  1,...,n  p,
where Y has multivariate joint normal distribution.
Is the set of feasible solutions convex when bounds are used?
Def: An m-multitree is called an m-multistar if the non-middle vertex
set of its multicherries are identical.
Rem: An m-multistar can be extended to an (m+1)-multistar that
provides us a better bound.
Th: Bounds based on multistars yield logarithmically concave
function in the probabilistic constraints, i.e.
h0* ( x)  BoundT1 x  Y1 ,...,Tn x  Yn 
is logarithmically concave if the correlations of Yj are
c1, j
cij 
ci ,i 1 for j  i  2,
c1,i 1
and c1 j ( j  2,...,n) and ci ,i 1 (i  2,...,n 1) are arbitrary positive numbers.
Exa: Covariance
cij 
i j for all 1 i < j.
REFERENCES
Bukszár, J. Upper Bounds for the Probability of a Union by Multitrees, Advances in Applied
Probability 33 (2), 437-452, 2001.
Bukszár, J. Prékopa, A. Probability Bounds with Cherry Trees, Mathematics of Operations
Research, 26 (1), 174-192, 2001.
Szántai, T. Bukszár, J. Probability Bounds given by Hypercherry Trees, Optimization Methods
and Software, 17 (3), 409-422, 2002.
Bukszár, J. Hypermultitrees and Bonferroni Inequalities, Mathematical Inequalities and
Applications, 6 (4), 727-743,2003.
Galambos, J. Simonelli, I. Bonferroni-type Inequalities with Applications, Springer-Verlag, NY, 1996.
Genz, A Numerical Computation of the Multivariate Normal Probabilities, J. Comput. Graph.
Stat. 1,141-150, 1992.
Grable, DA. Sharpened Bonferroni Inequalities, J. Combin. Theory Ser. B 57, 131-137, 1993.
Hoppe, FM., Seneta, E. A Bonferroni-type Identity and Permutation Bounds, International Statistical
Review 58, 3, 253-261, 1990.
Hunter, D. An Upper Bound for the Probability of a Union, J. Appl. Prob. 13, 597-603, 1976.
Prékopa, A. Stochastic Programming, Kluwer Academic Publishers, Dordrecht, 1995.
Prékopa, A. Boole-Bonferroni Inequalities and Linear Programming, Oper. Res. 36, 145-162,1988.
Prékopa, A. Sharp Bounds on Probabilities Using Linear Programming, Oper. Res. 38, 227-239,
1990.
Prékopa, A. The Discrete Moment Problem and Linear Programming, Discrete Applied
Mathematics, 27, 235-254, 1990.
Sobel, M. Uppuluri, VRR On Bonferroni-type Inequalities of the Same Degree for the Probability
of Unions and Intersections, Ann. Math. Statist. 43, 1549-1558, 1972.
Tomescu, I. Hypertrees and Bonferroni Inequalities, J. Combin. Theory Ser. B 41, 209-217, 1986.
Worsley, KJ. An Improved Bonferroni Inequality and Applications, Biometrika 69, 297-302, 1982.