import numpy as np
import numpy.random as rd
import scipy.stats as stats

def mypearson(x,y):
  mx=np.mean(x)
  my=np.mean(y)
  a=np.sum((x-mx)*(y-my))
  bx=np.sqrt(np.sum((x-mx)**2))
  by=np.sqrt(np.sum((y-my)**2))
  r=a/(bx*by)
  return r






mysize=1000000
data1=rd.poisson(50,mysize)
data2=rd.normal(loc=50,scale=10,size=mysize)

print "my pearson: ",mypearson(data1,data2)
print "numpy pearson: ",np.ma.corrcoef(data1,data2)[1,0]
print "scipy pearson: ",stats.pearsonr(data1,data2)[0]