Large Deviations of Bivariate Gaussian Extrema
Remco van der Hofstad Harsha Honnappa


Eindhoven University of Technology, P.O. Box 513, 5600 MB Eindhoven, The Netherlands. E-mail: rhofstad@win.tue.nl
Purdue University, 315 N. Grant St., West Lafayette IN 47906 U.S.A. Email: honnappa@purdue.edu
[
    Harsha Honnappa
    December 30, 2023
=====================





We establish sharp tail asymptotics for component-wise extreme values of bivariate Gaussian random vectors with arbitrary correlation between the components. We consider two scaling regimes for the tail event in which we demonstrate the existence of a restricted large deviations principle, and identify the unique rate function associated with these asymptotics. Our results identify when the maxima of both coordinates are typically attained by two different vs. the same index, and how this depends on the correlation between the coordinates of the bivariate Gaussian random vectors.

Our results complement a growing body of work on the extremes of Gaussian processes. The results are also relevant for steady-state performance and simulation analysis of networks of infinite server queues.





§ INTRODUCTION
Motivated by applications to the analysis of queueing networks, we study the large deviations of extreme values of multivariate Gaussians. We focus on the bivariate case for simplicity, but our analysis will carry over to the more general case with some effort. Let {X_1, X_2,…} be an ensemble of independent and identically distributed (i.i.d.) bi-variate Gaussians with covariance matrix Σ,
	and let X̅_n := (max_1≤ i≤ n X^(1)_i, max_1≤ i≤ nX_i^(2)) be the component-wise maximum, or extreme value, random vector. For simplicity we assume that [X_1] = 0. In the context of a queueing network X̅_n is an approximation to the maximum congestion experienced over a typical interval in a network of infinite server queues, for instance. We characterize the likelihood of the tail event {X̅_n > a_n u} where u ∈ (0,∞) as the typical interval n tends to infinity, under the assumption that a_n →∞ as n →∞. We consider two cases. 
	
Case 1: The right scale. Under the condition that a_n = √(log n), we prove a “restricted” large deviations principle (RLDP) (in the sense of <cit.>) in Theorem <ref> that shows that if u > √(2) (σ^(1),σ^(2)) (where σ^(j) is the standard deviation of marginal j) then
	
    lim_n→∞1/a_n^2log(X̅_n > a_n u) = J(u/σ),













where

    J(u)= 
     		1 - 1/2( u^(1))^2    when u^(2)≤ρ u^(1),
    
    		1 - 1/2( u^(2))^2    when u^(1)≤ρ u^(2),
    max{2 - 12u_2^2, 1-(u^(1)/σ^(1))^2-2ρ (u^(1)/σ^(1)) (u^(2)/σ^(2))+(u^(2)/σ^(2))^2/2(1-ρ^2)}   otherwise,

 u/σ := (u^(1)/σ^(1), u^(2)/σ^(2)), u_2^2 := ∑_i=1,2 |u^(i)|^2 and ρ∈ [-1,+1]. Here, the top two cases only arise when ρ>0, so when ρ<0, the last line equals J(u). The proof follows by using the Laplace principle in the key Lemma 1, and combining it with the `largest probability wins' principle. 

The different cases in (<ref>) originate due to the different scenarios in which the bivariate distribution can attain its maximum. In all the cases where a term +1 is present, the maximum is attained by one index of X_i which simultaneously attains the maximum of both coordinates. In all the cases where a term +2 is present, the maximum is attained by two different indices of X_i, one which attains the maximum of the first coordinate, and one which attains the maximum of the second coordinate. The latter case has most distinct possibilities (“larger entropy”), while the first may have a larger probability for appropriate correlation coefficients ρ. The optimal strategy is characterized by the `largest probability wins' principle. 
 	
Case 2: Larger scales. On a much larger scale, where a_n ≫√(log n), we establish two main results. First, we prove a leading order asymptote for the extreme value that aligns with the result in case 1. Precisely, in Theorem <ref> we prove the RLDP 
	
	
	lim_n→∞ 1/a_n^2 log(X̅_n> a_n u) = -I(u/σ),
	
where
	
	
	I(u)=
	1/2(u^(1))^2    when u^(2)≤ρu^(1),

	1/2 ( u^(2) )^2    when u^(1) ≤ρu^(2),

	1/2min{u^2_2, (u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2)/1-ρ^2}
	   otherwise.
	
	





In Theorem <ref>, on the other hand, we establish a sharp asymptote for the likelihood and show that there exists a continuous function I ^d → and constants K, b and c such that
	
    lim_n→∞ a_n^b n^c e^a_n^2 I(u)(X̅_n > a_n u) = K.

	The proofs of these theorems uses the inclusion-exclusion principle to bound the likelihood from above and below. 

	

  
Related Literature
Multivariate Gaussians emerge as stationary limits of networks of G/G/∞ infinite server queues when the arrival rate is high (i.e., in heavy-traffic in the sense of <cit.>); see <cit.> as well. This is straightforward to observe in the case of a single M/G/∞ station where the number in system in steady-state is Poisson distributed. When the arrival rate is high the steady-state distribution is well approximated by a Gaussian.


Next, there is an explicit connection with extreme value theory (EVT). The logarithmic asymptotics established here complement the uniform convergence results for EVT; see <cit.> and <cit.>.
There are also clear connections with recent work on extremes of multidimensional Gaussian processes in <cit.> and other related work, where logarithmic asymptotics are derived for the “at least one in the set" extremum (not the component-wise extremum considered here) for Gaussian processes. We note, in particular,  <cit.> where logarithmic asymptotics are derived for the “at least one in set" extremum of a sequence of (non-i.i.d.) generally distributed random vectors. The authors present a general theory closely aligned with the RLDP for univariate random variables introduced in <cit.>, whereby the Gärtner-Ellis condition need not be satisfied. Of course, our results are more restrictive in the sense that we only study i.i.d. Gaussian random vectors, but we also consider large-scale asymptotics that are not under consideration there.

Our results are also closely related to the important series of papers by Hashorva and Hüsler <cit.> generalizing the classic Mills ratio Gaussian tail bound <cit.>. We observe that the quadratic program logarithmic asymptote derived in Lemma 1 is also implied by the tail bound dervied in <cit.>. In <cit.>, the authors derive exact asymptotics for integrals of Gaussian random vectors, and in particular focus on the “at least one in the set" extremum for half-space extreme value sets. Our proof does not rely on the bound in <cit.>.

It would be interesting to strengthen Theorem <ref> to sharp asymptotics as is performed for large scales in Theorem <ref>. This is hard, since various error terms that can easily be dealt with in the proof of Theorem <ref> as they are much smaller than the leading order, will only become marginally smaller. In would also be interesting to extend our analysis to other multivariate random vectors with non-trivial dependence.



  
Notations and Setting
 All vector relations should be understood component-wise. Thus, x > y implies that x^ (j) > y^(j) for every component j. Following <cit.> we define a restricted large deviation principle (RLDP) as follows: for any q ∈^d a sequence of multivariate ^d-valued random variables {W_n} satisfies a RLDP with rate function J ^d → [0,∞] if 
	
    lim_n→∞1/v_nlog(W_n/a_n > q) = -J(q),

where v_n,a_n →∞ as n→∞. This asymptotic is not a full-fledged large deviations principle (LDP) since it does not provide any insight into what happens for negative q, i.e., it only deals with attaining large positive values. Furthermore, as noted in <cit.> and <cit.>, if W_n satisfies a LDP with  continuous rate function, then it automatically satisfies the RLDP. On the other hand, if the rate function is discontinuous, then it might not satisfy the RLDP.

We can write
	
	(X^(1)_i, X^(2)_i)d= (σ^(1)Z^(1)_i + μ^(1), σ^(2)Z^(2)_i + μ^(2)),
	
where (σ^(1),μ^(1)) and (σ^(2),μ^(2)) are the standard deviation and mean of X^(1)_i and  X^(2)_i, respectively, while (Z^(1)_i, Z^(2)_i) are standard bi-variate normals with correlation coefficient ρ∈[-1,1]. Assume that μ^(i) = 0 for i ∈{1,2}, without loss of generality.



§ RIGHT SCALE ASYMPTOTE

 We start by analyzing extreme events for bivariate Gaussian random variables:



Let {a_n}_n≥ 1 be any unbounded increasing sequence in n ∈. Then
  	
    lim_n→∞1/a_n^2log( X_1 > a_n ) = -_x > 1/2 x^T Σ^-1 x.

	

By definition, and with C an  explicit constant,
	
  	
  	1/a_n^2log( X_1 > a_n )    = 
	 1/a_n^2log(C ∫_x > a_n  exp(-1/2 x^TΣ^-1x) dx )

	    =1/a_n^2log( a_n C ∫_x >  exp(-a_n^21/2 x^TΣ^-1x) dx ),
where the second equality follows by a substitution of variables. Laplace's principle <cit.> implies the claim.


Next, we consider the asymptotics of the logarithmic likelihood of the event {∃ i ≤ n X_i > a_n u}. Note that this is not the component-wise maximum. 




Let a_n := √(log n). The bivariate Gaussian ensemble satisfies the RLDP limit
  	
    lim_n→∞1/a_n^2log(∃ i ≤ n  X_i > a_n u) = 1    - _x > u1/2 x^T Σ^-1 x,

for u := (u^(1),u^(2)) > √(2)(^(1),^(2)).




The condition (u^(1),u^(2)) > √(2)(^(1),^(2)) is very natural. Indeed, since the marginal distributions of each of the coordinates is normal with mean zero and standard deviation ^(j), we have that max_i=1^n X_i^(j)/√(logn)√(2)^(j) for j=1,2. Thus, when u^(j)≤√(2)^(j), it is natural to assume that this event does not contribute to the asymptotics in Proposition <ref>. In particular, when (u^(1),u^(2)) ≤√(2)(^(1),^(2)), the limit in (<ref>) equals zero, whereas the right-hand side is strictly positive.




Observe that

	
  	(∃i ≤n X_i > a_n u)    = ( ∪_i=1^n {X_i > a_n})

  	   = 1 - (∩_i=1^n { X_i > a_n u }^c)

  	   =(1-({X_1> a_n u}^c))∑_i=0^n-1 ({X_i > a_n u}^c)

  	   = (X_1 > a_n u) ∑_i=0^n-1 b_n^i,	 
	
where b_n := ({X_1 > a_n u}^c). 







From (<ref>), it follows that
  	
    log(∃ i ≤ n X_i > a_n u) ≤log(X_1 > a_n u) + log n,

using the fact that b_n < 1 for all finite n. Lemma <ref> implies that
  	
    lim sup_n→∞1/a_n^2log(X̅_n > a_n
       	 u) ≤ 1 - 1/2_x > u x^T Σ^-1 x.


Next, for the lower bound, we work with the term log∑_i=0^n-1 b_n^i to obtain a finer analysis. In particular,
suppose we demonstrate that, since b_n ∈ [0,1], log(n b_n^n-1) ≥log n + o(log n) as n →∞; then, it follows that
	
    log∑_i=0^n-1 b_n^i > log (n b_n^n-1) ≥log
      	n + o(log n) as n→∞.

Consequently, Lemma <ref>, combined with this result, implies that
	
    lim inf_n→∞1/a_n^2log( ∃ i ≤ n  X_i > a_n u) ≥
      	1 - _x > u1/2 x^T
      	Σ^-1 x,

thereby completing the proof of the proposition.

It remains to show (<ref>). Observe that the inclusion-exclusion formula implies that 
	
	
	b_n = ({X_1 > a_n u}^c)=({X_1^(1) ≤a_n u^(1)}∪{X_1^(2) ≤a_n u^(2)})
	
 satisfies
	
    b_n    = (X_1^(1)≤ a_n u^(1)) + (X_n^(2)≤ a_n u^(2)) - ({X_n^(1)≤ a_n u^(1)}∩{X_n^(2)≤ a_n u^(2)})
       = 2 - (X_1^(1) >  a_n u^(1)) - (X_n^(2) > a_n u^(2)) - ({X_n^(1)≤ a_n u^(1)}∩{X_n^(2)≤ a_n u^(2)}).

	
Therefore
	
	
	b_n ≥1 - (X_1^(1) >  a_n u^(1)) - (X_n^(2) > a_n u^(2)).
	
	
By the Taylor series expansion of log (1-x) = -x + o(x), as well as 1-b_n=o(1),
	
	
    log(n b_n^n-1)    = log n + (n-1) logb_n=log n - (n-1) (1-b_n) + o(log n)
       ≥log n +n(- (X_1^(1) >  a_n u^(1)) - (X_n^(2) > a_n u^(2)))(1 + o(1)).

	
Next, the Gaussian upper tail bound implies that for large n,
	
    log(n b_n^n-1)    ≥log n - n/√(2π)1/√(log n)∑_j∈{1,2}1/n^1/2(u^(i)/σ^(i))^21/σ^(j)u^(j)(1 + o(1))
       ≥log n - 1/√(2π)2/√(log n)max{1/n^(u^(1)/σ^(1))^2/2-11/σ^(1)u^(1), 1/n^(u^(2)/σ^(2))^2/2-11/σ^(2)u^(2)}(1 + o(1)).

Since (u^(1),u^(2)) > √(2)(^(1),^(2)), it follows that
	
	log(n b_n^n-1)/logn ≥1 + o(1)  as n→∞,
	
thereby completing the proof.

































	By a straightforward calculation,
	
	
	J_1(u/σ)   :=1-12_x > u x^T Σ^-1 x
	=
				1-1/2(u^(1)/σ^(1))^2    when u^(2)/σ^(2)≤ρu^(1)/σ^(1),

				1-1/2(u^(2)/σ^(2))^2    when u^(1)/σ^(1)≤ρu^(2)/σ^(2),

	1-(u^(1)/σ^(1))^2-2ρ(u^(1)/σ^(1)) (u^(2)/σ^(2))+(u^(2)/σ^(2))^2/2(1-ρ^2)
	   otherwise.
		





 Fix ρ∈ (0,1] and without loss of generality assume that σ^(1) = σ^(2) = 1. We can divide the positive quadrant into three regions as shown in Figure <ref>, where ρ = 0.5. Suppose that u is such that u^(2)≤ρ u^(1) (see the red region in Figure <ref>), then
	
    (u^(1)/σ^(1))^2-2ρ (u^(1)/σ^(1)) (u^(2)/σ^(2))+(u^(2)/σ^(2))^2/2(1-ρ^2)≤12(u^(1))^2 < 12(x^(1))^2,

where the final inequality holds for any x > u. It follows that 12_x > u x^T Σ^-1 x  = 12(u^(1))^2. A similar argument shows that 12_x > u x^T Σ^-1 x = 1/2(u^(2))^2 when u^(1)≤ρ u^(2). Finally, in the region where neither of these conditions holds (the blank region in Figure <ref>), it is straightforward to verify the Karush-Kuhn-Tucker (KKT) conditions for u and, since Σ^-1 is positive definite, u is the unique optimizer.



















As a consequence, we obtain the main result of this section:


Under the conditions of Proposition <ref>, 
	
	
	lim_n→∞ 1/a_n^2 log( X̅_n > a_n u) = J(u/σ),
	
where
	
    J(u)= 
     		1 - 1/2( u^(1))^2    when u^(2)≤ρ u^(1),
    
    		1 - 1/2( u^(2))^2    when u^(1)≤ρ u^(2),
    max{2 - 12u_2^2, 1-(u^(1)/σ^(1))^2-2ρ (u^(1)/σ^(1)) (u^(2)/σ^(2))+(u^(2)/σ^(2))^2/2(1-ρ^2)}   otherwise.

Further, with (I^⋆, J^⋆) the indices that maximize X̅_n (i.e., X̅_n=(X^(1)_I^⋆, X^(2)_J^⋆)),
	
	
	lim_n→∞(I^⋆≠J^⋆|X̅_n> a_n u)=
	
	1     when 2-12u/σ^2_2>J_1(u/σ),

	0     when 2-12u/σ^2_2<J_1(u/σ).
	
	



It is an interesting problem to extend (<ref>) to the case where 
	
	J(u/σ)=2-12u/σ^2_2=1-(u^(1)/σ^(1))^2-2ρ(u^(1)/σ^(1)) (u^(2)/σ^(2))+(u^(2)/σ^(2))^2/2(1-ρ^2),
	
but this seems quite difficult, as it requires a finer analysis beyond the principle of the largest term <cit.>.


Note that
	
	( X̅_n > a_n u)
	   =(∃i ≤nX_i > a_n u)

	     +(∃i≠j  ≤nX_i^(1) > a_n u^(1), X_j^(2) > a_n u^(1), ∄i ≤nX_i > a_n u).
Depending on u, the first or the second term will be dominant. Ignoring the event that ∄i ≤ n X_i > a_n u (which leads to an upper bound, but this event can also be incorporated in more detail), we have (using that X_i^(1), X_j^(2) are Gaussian)
	
	(∃i≠j  ≤nX_i^(1) > a_n u^(1), X_j^(2) > a_n u^(1))
	   ≈n(n-1) (X_i^(1) > a_n u^(1))(X_j^(2) > a_n u^(2))

	   ≈n^2 exp{-a_n^2 [(u^(1)/σ^(1))^2+(u^(2)/σ^(2))^2]/2}.
Taking log's and dividing by logn=a_n^2 gives 
	
	lim_n→∞ 1/a_n^2 log(∃i≠j  ≤nX_i^(1) > a_n u^(1), X_j^(2) > a_n u^(1))
	   =2-12[(u^(1)/σ^(1))^2+(u^(2)/σ^(2))^2]=2 - 12 u/σ_2^2

	   :=J_2(u/σ).
























Then, by the principle of the largest term <cit.>, it follows that
	
	
	lim_n→∞ 
	1/a_n^2 log(X̅_n > a_n u)    = max{lim_n→∞ 1/a_n^2 log( ∃i ≤nX_i > a_n u),

	           lim_n→∞ 1/a_n^2 log(∃i≠j  ≤nX_i^(1) > a_n u^(1), X_j^(2) > a_n u^(1) ) }

	   = max{J_1(u/σ), J_2(u/σ)}.
	
Further,
	
	lim_n→∞1/a_n^2 log(I^⋆≠J^⋆|X̅_n> a_n u)
	   =lim_n→∞1/a_n^2 [log(I^⋆≠J^⋆,X̅_n> a_n u)-log(X̅_n> a_n u)]

	   ≤J_2(u/σ)-J(u/σ)<0,
when J_2(u/σ)<J(u/σ), showing that (I^⋆≠ J^⋆|X̅_n> a_n u)=o(1) when J_2(u/σ)<J(u/σ).
Similarly,
	
	lim_n→∞1/a_n^2 log(I^⋆=J^⋆|X̅_n> a_n u)
	   =lim_n→∞1/a_n^2 [log(I^⋆=J^⋆,X̅_n> a_n u)-log(X̅_n> a_n u)]

	   ≤J_1(u/σ)-J(u/σ)<0,
when J_1(u/σ)<J(u/σ), showing that (I^⋆=J^⋆|X̅_n> a_n u)=o(1) when J_1(u/σ)<J(u/σ). This proves (<ref>) subject to (<ref>).


We are left to prove the claim in (<ref>), for which we consider several cases:

Case (1): u^(2)/σ^(2)≤ρ u^(1)/σ^(1). Under the assumption that u^ (j) > √(2)σ^ (j) for j = 1,2, it is straightforward to see that 1 - 1/2(u^(1)/σ^(1))^2 > 2 - 1/2u/σ_2^2. Lemma <ref> implies that
	
	
	lim_n→∞ 1/a_n^2 log(X̅_n > a_n u) = 1 - 12 (u^(1)/σ^(1))^2=J_1(u/σ).
	

Case (2): u^(2)/σ^(2)≤ρ u^(1)/σ^(1). The proof follows case (1) and is omitted.
	

Case (3): u^(2)/σ^(2) > ρ u^(1)/σ^(1) and u^(1)/σ^(1) > ρ u^(2)/σ^(2). Lemma <ref> implies that (<ref>) is simply
	
    max{J_1(u/σ), J_2(u/σ)}
    	=max{2 - 12u_2^2, 1-(u^(1)/σ^(1))^2-2ρ (u^(1)/σ^(1)) (u^(2)/σ^(2))+(u^(2)/σ^(2))^2/2(1-ρ^2)}.

Which of the two is the maximizer depends sensitively on the relation between ρ and u.

































































§ LARGE-SCALE ASYMPTOTE


  We look at (X̅_n>a_n u), where a_n≫√(logn), so that we are considering a large-deviation event. Recall that (X_i^(1),X^(2)) = (σ^(1) Z_i^(1), σ^(2) Z_i^(2)). By changing u=(u^(1), u^(2)) if needed, it thus suffices to study the standard case, and in what follows, we will therefore focus on (Z̅_n>a_n u). We prove the following two main theorems:



For any a_n≫√(logn), and with u^(2)≤ u^(1),
	
	
	 lim_n→∞ 1/a_n^2 log(Z̅_n> a_n u) = -I(u),
	
where 
	
	
	I(u)=
	1/2(u^(1))^2    when u^(2)≤ρu^(1),

	1/2min{u^2_2, (u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2/1-ρ^2}
	   otherwise.
	
	



Theorem <ref> is a special example of the two-dimensional Cramér Theorem (see e.g., <cit.> or <cit.>). The main aim of this paper in the large-scale asymptotics is to study the sharp asymptotics of a large deviation of the bi-variate normal distribution, which is quite interesting indeed:




For any a_n≫√(logn), and with u^(2)≤ u^(1),
	
	
	lim_n→∞ a_n^bn^c^a_n^2 I(u)(Z̅_n> a_n u) = K,
	
where 
	
	
	b=1+u^(2)≥ρu^(1),
	    c=1+I(u)=u^2/2,
	
and 
	
	K=
	
	1/2πu^(1)u^(2) 	    when I(u)=12u^2, u^(1)≠u^(2),

	1/4πu^(1)u^(2) 	    when I(u)=12u^2, u^(1) = u^(2),

	1/2πu^(1)				    when I(u)<12u^2, u^(2)<ρu^(1),

	1/4πu^(1)				    when I(u)<12u^2, u^(2)=ρu^(1),

	1-ρ^2/2π(u^(1)-ρu^(2))(u^(2)-ρu^(1))    when I(u)<12u^2, u^(2)>ρu^(1).
	
	

Consequently, with (I^⋆, J^⋆) the indices that maximize Z̅_n (i.e., Z̅_n=(Z^(1)_I^⋆, Z^(2)_J^⋆)),
	
	lim_n→∞(I^⋆≠J^⋆|Z̅_n> a_n u)=
	
	1     when I(u)=12u^2_2,

	0     otherwise.
	
	



We prove Theorem <ref> and <ref> in one go. We note that
	
	(Z̅_n>a_n u)=(⋃_(i,j) {Z^(1)_i>a_n u^(1), Z^(2)_j> a_n u^(2)}).
	
We obtain
	
	(Z̅_n>a_n u)=(⋃_(i,j) A_(i,j)),
	
where 
	
	A_(i,j)={Z^(1)_i>a_n u^(1), Z^(2)_j> a_n u^(2)}.
	
From this formula, we see the importance of symmetry, as A_(i,j)=A_(j,i) for the symmetric case where u^(1)=u^(2), but not when this is not the case. 
In the symmetric case where u^(1)=u^(2), we note that A_(i,j)=A_(j,i), so we may write, instead,
	
	(Z̅_n>a_n u)=(⋃_(i,j)i≤j {Z^(1)_i>a_n u^(1), Z^(2)_j> a_n u^(2)}).
	



  
Using inclusion-exclusion. We use inclusion-exclusion to obtain that
	
	∑_(i,j) (A_(i,j))
	-e_n(u)≤(Z̅_n>a_n u)≤∑_(i,j) (A_(i,j)),
	
where
	
	
	e_n(u)=1/2∑_(i,j)≠(k,l)(A_(i,j)∩A_(k,l)),
	
while for the symmetric case, we sum over ordered pairs (i,j) with i≤ j instead.

Below, we analyse each of these terms. We separate between the case where (i) the indices are different; (ii) they are equal but the probability simplifies; (iii) they are different and we need to perform the integral over the joint density using the Laplace method. We start with the asymmetric case where u^(1)≠ u^(2), remarking on the extension to the symmetric case at the end of the proof. Without loss of generality, we may assume that u^(1)>u^(2).





  
Sum of probabilities: unequal indices.  First consider the case where i≠ j. Then, since (Z^(1)_i, Z^(2)_j) are i.i.d. standard normal random variables, 
	
	∑_(i,j)i≠j (A_(i,j))
	=n(n-1) [1-Φ(a_n u^(1))][1-Φ(a_n u^(2))],
	
where Φ(x)=(Z≤ x) is the error-function or the distribution function of a standard normal. By the asymptotics, for x large,
	
	1-Φ(x)=1/√(2π) x^-x^2/2(1+O(x^-2)),
	
we thus obtain that
	
	
	∑_(i,j)i≠j (A_(i,j))
	=n^2 a_n^2 1/2πu^(1)u^(2)^-a_n^2 u_2^2/2(1+o(1)).
	



  
Sum of probabilities: simple cases of equal indices. We next proceed with the case where i=j, for which we get	
	
	∑_(i,i) (A_(i,i))
	=n(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2)).
	
We use that, conditionally on Z^(1), the law of Z^(2) equals Z^(2)=ρ Z^(1) +√(1-ρ^2)Z, where Z is independent of Z^(1). We thus get that
	
	(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))
	   =(Z^(1)>a_n u^(1), ρZ^(1) +√(1-ρ^2)Z> a_n u^(2))

	   =[Z^(1)>a_n u^(1)(√(1-ρ^2) Z>a_n u^(2)-ρZ^(1)|Z^(1) )].
	
When ρ u^(1)>u^(2),
	
	(√(1-ρ^2) Z>a_n u^(2)-ρZ^(1)|Z^(1) )=1+o(1),
	
so that 
	
	
	∑_(i,i) (A_(i,i))
	=n(Z^(1)>a_n u^(1))(1+o(1))
	=na_n 1/2πu^(1)^-a_n^2 (u^(1))^2/2(1+o(1)).
	
When ρ u^(1)=u^(2), instead
	
	(√(1-ρ^2) Z>a_n u^(2)-ρZ^(1)|Z^(1) )=1/2+o(1),
	
since Z^(1)-a_n u^(1)=(1) when Z^(1)>a_n u^(1). Thus, for ρ u^(1)=u^(2) this leads to 
	
	
	∑_(i,i) (A_(i,i))
	=n(Z^(1)>a_n u^(1))(1+o(1))
	=na_n 1/4πu^(1)^-a_n^2 (u^(1))^2/2(1+o(1)).
	



  
Sum of probabilities: Lapace integral for equal indices. When the above simple cases do not apply, we write (Z^(1)>a_n u^(1), Z^(2)> a_n u^(2)) explicitly as a two-dimensional integral as
	
	(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))
	=1/2π∫_a_n u^(1)^∞ ∫_a_n u^(2)^∞  ^-(x_1^2 -2ρx_1x_2 +x_2^2)/2(1-ρ^2)dx_2dx_1.
	
We rescale the integrands by a_n to obtain
	
	(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))
	=a_n^2/2π∫_u^(1)^∞ ∫_u^(2)^∞  ^-a_n^2(x_1^2 -2ρx_1x_2 +x_2^2)/2(1-ρ^2)dx_2dx_1.
	
This is a classic example of a Laplace integral. Thus, the integral is dominated by the minimum of (x_1^2 -2ρ x_1x_2 +x_2^2)/2(1-ρ^2) over all (x_1,x_2) for which x_1≥ u^(1), x_2≥ u^(2). Since (x_1^2 -2ρ x_1x_2 +x_2^2)/2(1-ρ^2) is convex, this minimum is attained at one of the boundaries. Since ρ u^(1)<u^(2), this minimum is attained at x_1=u^(1), x_2=u^(2) (see also the analysis in Lemma <ref>). 


Thus,
	
	(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))
	   =a_n^2/2πexp{-a_n^2 (u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2/2(1-ρ^2)}

	   ×∫_u^(1)^∞ ∫_u^(2)^∞  exp{-a_n^2(x_1^2 -2ρx_1x_2 +x_2^2)-(u^(1))^2+2ρu^(1) u^(2)-(u^(2))^2/2(1-ρ^2)dx_2dx_1.
	
Therefore, we obtain that
	
	   2πa_n^-2 exp{a_n^2 (u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2/2(1-ρ^2)}(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))

	   =∫_u^(1)^∞ ∫_u^(2)^∞  exp{-a_n^2(x_1-u^(1))^2 -2ρ(x_1-u^(1))(x_2-u^(2)) +(x_2-u^(2))^2/2(1-ρ^2)

	     ×exp{-a_n^2[2u^(1)(x_1-u^(1))+2u^(2)(x_2-u^(2))
	-ρu^(1)(x_2-u^(2))-ρu^(2)(x_1-u^(1))]/2(1-ρ^2)}dx_2dx_1.
	
Since ρ u^(1)<u^(2), we have that the quadratic function inside the exponential is minimized for x_1=u^(1), x_2=u^(2). Shifting both integrands by u^(1) and u^(2) respectively, leads to
	
	   2πa_n^-2 exp{a_n^2 (u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2/2(1-ρ^2)}(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))

	   =∫_0^∞ ∫_0^∞  exp{-a_n^2x_1^2 -2ρx_1x_2 +x_2^2/2(1-ρ^2)}exp{-a_n^22u^(1)x_1+2u^(2)x_2
	-ρu^(1)x_2-ρu^(2)x_1/2(1-ρ^2)}dx_2dx_1.
	
Now rescaling both integrands by a_n^-2 leads to 
	
	   2πa_n^2 exp{a_n^2 (u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2/2(1-ρ^2)}(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))

	   =∫_0^∞ ∫_0^∞  exp{-a_n^-2x_1^2 -2ρx_1x_2 +x_2^2/2(1-ρ^2)}
	exp{-(u^(1)-ρu^(2))x_1+(u^(2)-ρu^(1))x_2/1-ρ^2}dx_2dx_1.
	
Again we see the significance of the assumption that ρ u^(1)<u^(2), which implies that both linear terms have a negative coefficient, and thus the exponential functions are integrable. As a result, the first exponential in the integral only leads to an error term, so that 
	
	
	   2πa_n^2 exp{a_n^2 (u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2/2(1-ρ^2)}(Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))

	   =(1+o(1))∫_0^∞ ∫_0^∞exp{-(u^(1)-ρu^(2))x_1+(u^(2)-ρu^(1))x_2/1-ρ^2}dx_2dx_1

	   =(1+o(1)) (1-ρ^2)(u^(1)-ρu^(2))^-1(u^(2)-ρu^(1))^-1.
	
Combining (<ref>)–(<ref>) with (<ref>) yields the asymptotics of the sum of probabilities. Note that the final outcome yields (<ref>) in the asymmetric case, so what is left is to show that the error term e_n(u) is of smaller order.




  
The symmetric case: sum of probabilities. We now look at the sum of probabilities for the symmetric case, and analyse (A_(i,j)) there. The analysis for the case where i≠ j is identical to the one above, except for the fact that the prefactor (due to the number of pairs (i,j)) is changed from n(n-1) to n(n-1)/2. The contribution for the case where i=j is also the same as above. In fact, it is easy to see that for u^(1)=u^(2)=u, we have I((u,u))=u^2 when ρ≤ 0, while I((u,u))=u^2/(1+ρ) for ρ>0, since
	
	(u^(1))^2-2ρu^(1) u^(2)+(u^(2))^2/1-ρ^2
	=u^2 2(1-ρ)/1-ρ^2=u^2 2/1+ρ<2u^2,
	
precisely when ρ>0. 



  
The error term e_n(u): asymmetric case.  In dealing with error terms, we will make essential use of the fact that a_n≫√(logn). This condition implies that if a certain event A satisfies (A)≤^-a_n^2 J for some J>I(u), then (A) will constitute an error term in evaluating (Z̅_n>a_n u), irrespective of the precise powers of a_n and n. Recall (<ref>). We investigate the different ways that (i,j)≠ (k,l) can occur, depending on the cardinality of {i,j,k,l} which ranges from 2 to 4. Below, we assume throughout the analysis that the indices i,j,k,l used are distinct.




  Case (2a): (i,j), (j,i) This corresponds to (Z̅_n>a_n u̅), where u̅=(u^(1)∨ u^(2), u^(1)∨ u^(2)) and x∨ y=max{x,y} for x,y∈ℝ. This case was investigated in the previous step, and we see that the rate at speed a_n^2 equals I((u̅,u̅))>I(u), since u^(1)≠ u^(2).


  Case (2b): (i,i), (i,j) or (i,i), (j,i) By independence, these probabilities equal (Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))(Z>a_nu^(2)) and  (Z^(1)>a_n u^(1), Z^(2)> a_n u^(2))(Z>a_nu^(1)), respectively. Obviously, the rate at speed a_n^2 is strictly larger than I(u).


  Case (3a): (i,i), (j,k) By independence, this probability equals (A_(i,i))(A_(j,k)), the rate at speed a_n^2 again being strictly larger than I(u).


  Case (3b): (i,j), (j,k) By independence, this probability equals (A_(j,j))(Z>a_nu^(1))(Z>a_nu^(2)), the rate at speed a_n^2 again being strictly larger than I(u).


  Case (3c): (i,j), (k,j) This case is similar. 


  Case (4): (i,j), (k,l) By independence, this probability equals (A_(i,j))(A_(k,l)), the rate at speed a_n^2 being at least 2I(u), which is again strictly larger than I(u).

Together, these cases show that e_n(u) is of smaller order than the sum of probabilities in the asymmetric case.



  
The error term e_n(u): symmetric case. 	This analysis is similar the the asymmmetric case, except that some cases do not arise. We again go through the distinct possibilities, writing u=u^(1)=u^(2):




  Case (2): (i,i), (i,j) or (i,i), (j,i) By independence, this probability equals (Z^(1)>a_n u, Z^(2)> a_n u)(Z>a_nu). Obviously, the rate at speed a_n^2 is strictly larger than I((u,u)).


  Case (3a): (i,i), (j,k) By independence, this probability equals (A_(i,i))(A_(j,k)), the rate at speed a_n^2 again being strictly larger than I((u,u)).


  Case (3b): (i,j), (j,k) By independence, this probability equals (A_(j,j))(Z>a_nu)(Z>a_nu), the rate at speed a_n^2 again being strictly larger than I((u,u)).


  Case (4): (i,j), (k,l) By independence, this probability equals (A_(i,j))(A_(k,l)), the rate at speed a_n^2 being at least 2I((u,u)), which is again strictly larger than I((u,u)).

Together, these cases show that e_n(u) is also of smaller order than the sum of probabilities in the symmetric case.










  
Acknowledgement. The work of RvdH is supported by the Netherlands Organisation for Scientific Research (NWO) through VICI grant 639.033.806 and the Gravitation Networks grant 024.002.003. The work of HH is partially supported by the National Science Foundation through grants CMMI-1636069 and DMS-1812197.

plain



