docs: def calc_LD(M,pos,i_start=[0],max_dist=1000000):
      """calculate linkage disequilibrium correlations:
      M: genotype matrix
      pos: position vector
      i_start: index to start from for LD calculation
      dist: distance
      """
      RV = []
      DIST = []
      for start in i_start:
          pos0 = pos[start]
          v0  = M[:,start]
          Iselect = sp.nonzero(sp.absolute(pos-pos0)<=max_dist)[0]
          rv = sp.zeros(len(Iselect))
          for i in range(len(Iselect)):
              rv[i] = (sp.corrcoef(v0,M[:,Iselect[i]])[0,1])**2
          #sort by distance
          dist = sp.absolute(pos[Iselect]-pos0)
          RV.extend(rv)
          DIST.extend(dist)
      RV = sp.array(RV)
      DIST = sp.array(DIST)
      II = DIST.argsort()
      return [DIST[II],RV[II]]