Advertisement
Guest User

Untitled

a guest
Jun 19th, 2012
493
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 0.78 KB | None | 0 0
  1. skew normal distribution in scipy
  2. from scipy import linspace
  3. from scipy import pi,sqrt,exp
  4. from scipy.special import erf
  5.  
  6. from pylab import plot,show
  7.  
  8. def pdf(x):
  9. return 1/sqrt(2*pi) * exp(-x**2/2)
  10.  
  11. def cdf(x):
  12. return (1 + erf(x/sqrt(2))) / 2
  13.  
  14. def skew(x,e=0,w=1,a=0):
  15. t = (x-e) / w
  16. return 2 / w * pdf(t) * cdf(a*t)
  17. # You can of course use the scipy.stats.norm versions
  18. # return 2 * norm.pdf(t) * norm.cdf(a*t)
  19.  
  20.  
  21. n = 2**10
  22.  
  23. e = 1.0 # location
  24. w = 2.0 # scale
  25.  
  26. x = linspace(-10,10,n)
  27.  
  28. for a in range(-3,4):
  29. p = skew(x,e,w,a)
  30. plot(x,p)
  31.  
  32. show()
  33.  
  34. fzz = skew(x,e,w,a) + norm.rvs(0,0.04,size=n) # fuzzy data
  35.  
  36. def optm(l,x):
  37. return skew(x,l[0],l[1],l[2]) - fzz
  38.  
  39. print leastsq(optm,[0.5,0.5,0.5],(x,))
  40.  
  41. (array([ 1.05206154, 1.96929465, 0.94590444]), 1)
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement