importsympyassymimportnumpyasnpdeftheoretic_fixation(N,game,i=1):""" Calculate x_i as given by the above formula """f_ones=np.array([(game[0,0]*(i-1)+game[0,1]*(N-i))/(N-1)foriinrange(1,N)])f_twos=np.array([(game[1,0]*i+game[1,1]*(N-i-1))/(N-1)foriinrange(1,N)])gammas=f_twos/f_onesreturn(1+np.sum(np.cumprod(gammas[:i-1])))/(1+np.sum(np.cumprod(gammas)))