def gelman_rubin(mtrace):
 if (mtrace.nchains < 2):
  raise ValueError('Gelman-Rubin diagnostic requires multiple chains of the same length.')
 Rhat = {}
 for var in mtrace.varnames:
  x = np.array(mtrace.get_values(var, combine=False))
  num_samples = x.shape[1]
  B = (num_samples * np.var(np.mean(x, axis=1), axis=0, ddof=1))
  W = np.mean(np.var(x, axis=1, ddof=1), axis=0)
  Vhat = (((W * (num_samples - 1)) / num_samples) + (B / num_samples))
  Rhat[var] = np.sqrt((Vhat / W))
 return Rhat

def gelman_rubin ( mtrace ) : if ( mtrace . nchains < 2 ) : raise ValueError ( 'Gelman-Rubin diagnostic requires multiple chains of the same length.' ) Rhat = { } for var in mtrace . varnames : x = np . array ( mtrace . get_values ( var , combine = False ) ) num_samples = x . shape [ 1 ] B = ( num_samples * np . var ( np . mean ( x , axis = 1 ) , axis = 0 , ddof = 1 ) ) W = np . mean ( np . var ( x , axis = 1 , ddof = 1 ) , axis = 0 ) Vhat = ( ( ( W * ( num_samples - 1 ) ) / num_samples ) + ( B / num_samples ) ) Rhat [ var ] = np . sqrt ( ( Vhat / W ) ) return Rhat
returns estimate of r for a set of traces .