of ¯ ζ$_{t}$ is, σ $^{2}$[ ¯ ζ$_{t}$ ] = ∑ k = t 1 N$_{T}$$_{k}$ ( 〈 α 2 $_{k}$π ( x $^{′}$) $^{2}$〉$_{π}$$_{k}$$_{,T}$$_{k}$ − 〈 α$_{k}$π ( x $^{′}$) 〉 2 π$_{k}$,T$_{k}$ ) ∑ k = t 〈 α$_{k}$π ( x $^{′}$) 〉 2 π$_{k}$,T$_{k}$ + ∑ k = t 1 N$_{π}$$_{k}$ ( 〈 α 2 $_{k}$T$_{k}$ ( x $^{′}$| x ) $^{2}$〉$_{π}$$_{k}$$_{,π}$ − 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉 2 π$_{k}$,π$_{t}$ ) ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉 2 π$_{k}$,π$_{t}$ (9) = ∑ k = t 〈 α 2 $_{k}$T$_{k}$ ( x $^{′}$| x ) ( 1 $_{N$_{T}$$_{k}$}$π ( x $^{′}$) e − ζ + 1 $_{N$_{π}$$_{k}$}$T$_{k}$ ( x $^{′}$| x ) )〉 π$_{k}$,π$_{t}$ [ ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉$_{π}$$_{k}$$_{,π}$$_{t}$ ] 2 − ∑ k = t ( 1 N$_{T}$$_{k}$ + 1 N$_{π}$$_{k}$ ) 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉 2 π$_{k}$,π$_{t}$ [ ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉$_{π}$$_{k}$$_{,π}$$_{t}$ ] 2 (10) ∑ k = t ( 1 N$_{T}$$_{k}$ + 1 N$_{π}$$_{k}$ ) 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉 2 π$_{k}$,π$_{t}$ [ ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉$_{π}$$_{k}$$_{,π}$$_{t}$ ] 2 ≤ max ( 1 N$_{T}$$_{k}$ + 1 N$_{π}$$_{k}$ ) ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉 2 π$_{k}$,π$_{t}$ [ ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉$_{π}$$_{k}$$_{,π}$$_{t}$ ] 2 ≤ max ( 1 N$_{T}$$_{k}$ + 1 N$_{π}$$_{k}$ ) [ ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉$_{π}$$_{k}$$_{,π}$$_{t}$ ] 2 [ ∑ k = t 〈 α$_{k}$T$_{k}$ ( x $^{′}$| x ) 〉$_{π}$$_{k}$$_{,π}$$_{t}$ ] 2 = max ( 1 N$_{T}$$_{k}$ + 1 N$_{π}$$_{k}$ ) (11) ∑ k = t  $^{N$_{π}$$_{k}$}$ ∑ l =1 1 1 + N$_{π}$$_{k}$ N$_{T}$$_{k}$ π ( x ′ $_{l}$) T$_{k}$ ( x ′ $_{l}$| x$_{l}$ $_{)}$e − ¯ ζ − N$_{T}$$_{k}$ ∑ m =1 1 1 + N$_{T}$$_{k}$ N$_{π}$$_{k}$ T$_{k}$ ( x ′ $_{m}$| x$_{m}$ ) π ( x ′ $_{m}$) e ¯ $_{ζ}$  = 0 , (12) Unlike in the case where K = 2, the second term in the above variance expression does not reduce to a constant. We can show, however, using a procedure similar to Guibas and Veach [6, 7], that it is bounded by a constant: Thus, if we find a set of functions for α$_{k}$ that minimize the first term in σ [ ¯ ζ$_{t}$ ], no other choices for α$_{k}$ can decrease the variance by more than max ( 1 N$_{T}$$_{k}$ + 1 N$_{π}$$_{k}$ ) . Minimization of the first term in Eq. 9 by Lagrange multipliers leads to an optimized set of α$_{k}$ analogous to Eq. 4. As with K = 2, substitution of this α$_{k}$ leads to an estimator for the integral which can either be solved self-consistently or by finding the zero of an implicit function, The corresponding variance may be estimated by substituting the appropriate α$_{k}$ into Eq. 9. The original and optimized forms of the replica gas identity were demonstrated by esti- mating the partition function of a two-dimensional 32 × 32 spin Ising model, using parame- ters similar to those previously described [3]. As in a typical Ising model simulation without 5