src/daComposant/daDiagnostics/GaussianAdequation.py

   1 #-*-coding:iso-8859-1-*-
   2 #
   3 #  Copyright (C) 2008-2010  EDF R&D
   4 #
   5 #  This library is free software; you can redistribute it and/or
   6 #  modify it under the terms of the GNU Lesser General Public
   7 #  License as published by the Free Software Foundation; either
   8 #  version 2.1 of the License.
   9 #
  10 #  This library is distributed in the hope that it will be useful,
  11 #  but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 #  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 #  Lesser General Public License for more details.
  14 #
  15 #  You should have received a copy of the GNU Lesser General Public
  16 #  License along with this library; if not, write to the Free Software
  17 #  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  18 #
  19 #  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  20 #
  21 __doc__ = """
  22     Diagnostic qui effectue le test du Khi2 pour juger de l'adéquation entre
  23     la distribution d'un échantillon et une distribution gaussienne dont la
  24     moyenne et l'écart-type sont calculés sur l'échantillon.
  25     En input : la tolerance(tolerance) et le nombre de classes(nbclasse)
  26     En output : Le resultat du diagnostic est une reponse booleenne au test :
  27                 True si l adequation a une distribution gaussienne est valide
  28                 au sens du test du Khi2,
  29                 False dans le cas contraire.
  30 """
  31 __author__ = "Sophie RICCI - Juillet 2008"
  32
  33 import numpy
  34 from daCore import BasicObjects, Persistence
  35 from ComputeKhi2 import ComputeKhi2_Gauss
  36 import logging
  37
  38 # ==============================================================================
  39 class ElementaryDiagnostic(BasicObjects.Diagnostic,Persistence.OneScalar):
  40     """
  41     """
  42     def __init__(self, name="", unit="", basetype = None, parameters = {} ):
  43         BasicObjects.Diagnostic.__init__(self, name, parameters)
  44         Persistence.OneScalar.__init__( self, name, unit, basetype = bool)
  45         for key in ["tolerance", "dxclasse", "nbclasses"]:
  46             if not self.parameters.has_key(key):
  47                 raise ValueError("A parameter named \"%s\" is required."%key)
  48
  49     def formula(self, V):
  50         """
  51         Effectue le calcul de la p-value pour un vecteur et une distribution
  52         gaussienne et un nombre de classes donne en parametre du diagnostic.
  53         """
  54
  55         [vectclasse, eftho, efobs, valeurKhi2, areaKhi2, message] = ComputeKhi2_Gauss(
  56             vectorV = V,
  57             dx = self.parameters["dxclasse"],
  58             nbclasses = self.parameters["nbclasses"],
  59             SuppressEmptyClasses = True)
  60
  61         logging.info( message )
  62         logging.info( "(si <%.2f %s on refuse effectivement l'adéquation)"%(100.*self.parameters["tolerance"],"%") )
  63         logging.info("vecteur des classes=%s"%numpy.size(vectclasse) )
  64         logging.info("valeurKhi2=%s"%valeurKhi2)
  65         logging.info("areaKhi2=%s"%areaKhi2)
  66         logging.info("tolerance=%s"%self.parameters["tolerance"])
  67
  68         if (areaKhi2 < (100.*self.parameters["tolerance"])) :
  69             answerKhisquareTest = False
  70         else:
  71             answerKhisquareTest = True
  72         logging.info( "La réponse au test est donc est %s"%answerKhisquareTest )
  73         return answerKhisquareTest
  74
  75     def calculate(self, vector = None,  step = None):
  76         """
  77         Active la formule de calcul
  78         """
  79         if vector is None:
  80             raise ValueError("One vector must be given to calculate the Khi2 test")
  81         V = numpy.array(vector)
  82         if V.size < 1:
  83             raise ValueError("The given vector must not be empty")
  84         #
  85         value = self.formula( V )
  86         #
  87         self.store( value = value,  step = step)
  88
  89 # ==============================================================================
  90 if __name__ == "__main__":
  91     print "\n AUTODIAGNOSTIC \n"
  92
  93     print " Test d adequation du khi-2 a une gaussienne pour un vecteur x"
  94     print " connu de taille 1000, issu d'une distribution gaussienne normale"
  95     print " en fixant la largeur des classes"
  96     print
  97     #
  98     # Initialisation des inputs et appel du diagnostic
  99     # ------------------------------------------------
 100     tolerance = 0.05
 101     dxclasse = 0.1
 102     D = ElementaryDiagnostic("AdequationGaussKhi2", parameters = {
 103         "tolerance":tolerance,
 104         "dxclasse":dxclasse,
 105         "nbclasses":None,
 106         })
 107     #
 108     # Tirage de l'echantillon aleatoire
 109     # ---------------------------------
 110     numpy.random.seed(2490)
 111     x = numpy.random.normal(50.,1.5,1000)
 112     #
 113     # Calcul
 114     # ------
 115     D.calculate(x)
 116     #
 117     if D.valueserie(0) :
 118         print " L'adequation a une distribution gaussienne est valide."
 119         print
 120     else :
 121         raise ValueError("L'adéquation a une distribution gaussienne n'est pas valide.")
 122
 123
 124     print " Test d adequation du khi-2 a une gaussienne pour u:n vecteur x"
 125     print " connu de taille 1000, issu d'une distribution gaussienne normale"
 126     print " en fixant le nombre de classes"
 127     print
 128     #
 129     # Initialisation des inputs et appel du diagnostic
 130     # ------------------------------------------------
 131     tolerance = 0.05
 132     nbclasses = 70.
 133     D = ElementaryDiagnostic("AdequationGaussKhi2", parameters = {
 134         "tolerance":tolerance,
 135         "dxclasse":None,
 136         "nbclasses":nbclasses
 137         })
 138     #
 139     # Tirage de l'echantillon aleatoire
 140     # ---------------------------------
 141     numpy.random.seed(2490)
 142     x = numpy.random.normal(50.,1.5,1000)
 143     #
 144     # Calcul
 145     # ------
 146     D.calculate(x)
 147     #
 148     if D.valueserie(0) :
 149         print " L'adequation a une distribution gaussienne est valide."
 150         print
 151     else :
 152         raise ValueError("L'adequation a une distribution gaussienne n'est pas valide.")
 153
 154