Skip to Content

stats_circulaires.py

  1. import numpy as np
  2.  
  3. # solution 1 avec les cosinus directeurs
  4. # --------------------------------------
  5. from collections import namedtuple
  6. def cosdir_azim(azim):
  7. az = np.radians(azim)
  8. cosa =np.sin(az)
  9. cosb = np.cos(az)
  10. name = namedtuple("name", ["cosa", "cosb"])
  11. return name(cosa,cosb)
  12.  
  13. angles = np.array([341.0, 359.0, 334.0, 15.0, 330.0, 301.0, 299.0, 9.0, 7.0, 353.0, 28.0, 25.0, 23.0, 30.0, 350.0, 25.0, 22.0, 8.0, 356.0, 27.0])
  14. sommecosa = sum(cosdir_azim(angles).cosa)
  15. sommecosb = sum(cosdir_azim(angles).cosb)
  16.  
  17. # direction moyenne angulaire
  18. print "{:.6f}".format(np.degrees(np.arctan2(sommecosa,sommecosb)))
  19. 1.060902
  20. # longueur résultante moyenne
  21. print "{:.6f}".format(np.sqrt(sommecosa*sommecosa + sommecosb*sommecosb)/len(angles))
  22. 0.896237
  23.  
  24. # variance angulaire
  25. print 1 - np.sqrt(sommecosa*sommecosa + sommecosb*sommecosb)/len(angles)
  26. 0.103763188613
  27.  
  28. # ou directement
  29.  
  30. def getCircularMean(angles):
  31. n = len(angles)
  32. sineMean = np.divide(np.sum(np.sin(np.radians(angles))), n)
  33. cosineMean = np.divide(np.sum(np.cos(np.radians(angles))), n)
  34. vectorMean = np.arctan2(sineMean, cosineMean)
  35. return np.degrees(vectorMean)
  36.  
  37. print "{:.6f}".format(getCircularMean(angles))
  38. 1.060902
  39.  
  40. # solution 2 avec SciPy
  41. # ---------------------
  42.  
  43. from scipy import stats
  44. T1 = np.radians(np.array(angles))
  45. # direction moyenne angulaire
  46. print "{:.6f}".format(np.degrees(stats.circmean(T1)))
  47. 1.060902
  48.  
  49.  
  50. # solution 3 avec pycircstat
  51. # ---------------------------
  52.  
  53. import pycircstat
  54. T1 = np.radians(np.array(angles))
  55. # direction moyenne angulaire
  56. print "{:.6f}".format(np.degrees(np.array(pycircstat.mean(T1))))
  57. 1.060902
  58. # longueur résultante moyenne
  59. print "{:.6f}".format(pycircstat.resultant_vector_length(T1))
  60. 0.896237
  61. # variance angulaire
  62. 1 - pycircstat.resultant_vector_length(T1
  63. 0.1037632
  64. # variance angulaire selon Batschelet (1981) = 2 * (1 - R)
  65. print "{:.6f}".format(pycircstat.avar(T1))
  66. 0.207526
  67. # déviation angulaire standard
  68. print "{:.6f}".format(pycircstat.std(T1))
  69. 0.468082