MATH_IntegraleGamma    [Statistiques]

Par Teddy Linet, Philip Burns (Math4D v2)
Nouvelle recherche
Si (Faux)
   ` MATH_IntegralGammaIncomp [Teddy LINET 02/2003 et Philip BURNS 1985]
   ` ---------------------------------------------
   ` ATTENTION : nécessite MATH_LogGamma
   ` ---------------------------------------------
   ` Rappels mathématiques :
   ` Evalue l'intégrale de la fonction Gamma "Incomplète"
   ` soit de 0 à la valeur $y_F
   ` 
   ` La fonction gamma est l'intégrale sur t entre 0 et + l'infini de la fonction
   ` f(x)=(1/gamma(alpha))*(lambda^alpha)*(x^(alpha-1))*exp(-lambda*x)
   ` ou gamma(t) est l'intégrale entre 0 et + l'infini de f(t)=t^z-1*exp(-t)
   ` Elle a entre autre la particularité que gamma(z+1)=z*gamma(z)
   `  
   ` Elle est utile en statistiques car par exemple on peut démontrer que la
   ` distribution d'une fonction khi 2 suit une fonction gamma(X2/2;ddl/2)
   `
   ` Remarque : Utilise en fonction des valeurs des arguments
   ` soit une somme infinie de série soit une fraction continue
   ` ---------------------------------------------
   ` Ecriture :
   ` MATH_IntegralGammaIncomp($y_F; $ddl_F;$vDprec_I ;$vMaxIter_I)
   ` $y_F : Valeur de la distribution gamma
   ` $ddl_F : degrés de liberté 
   ` $vDprec_I : degré de précision
   ` $vMaxIter_I : Nombre max d'iteration
   ` ---------------------------------------------  
   ` MATHERROR
   ` 0 -> Pas d'erreur
   ` -1 -> $y_F <=0
   ` -2 -> $ddl_F <=0
   ` -3 -> La valeur du paramètre intermédiaire $vF_F est inférieure au minimum
   ` -4 -> La valeur du paramètre intermédiaire exp($vF_F) est nulle
   ` ---------------------------------------------
Fin de si 

C_REEL($1;$y_F;$2;$ddl_F;$vIntegrale_F)
C_ENTIER($3;$vDprec_I;$4;$vMaxIter_I)
C_REEL($vCprec_F)
C_ENTIER($vIter_I)
C_ENTIER LONG(MATHERROR)
$y_F:=$1
$ddl_F:=$2
  ` S'adapter en fonction du nombre de parametres
Si (Nombre de parametres>2)
 $vDprec_I:=Abs($3)
Sinon 
 $vDprec_I:=12
Fin de si 
Si (Nombre de parametres>3)
 $vMaxIter_I:=Abs($4)
Sinon 
 $vMaxIter_I:=200
Fin de si 
  `VARIABLES
C_REEL($vF_F;$vC_F;$vA_F)
C_REEL($B_F;$vTerm_F)
TABLEAU REEL($tPn_F;6)
C_REEL($vGin_F;$vAn_F;$vRn_F)
C_REEL($vDif_F;$vEps_F)
C_BOOLEEN($vDone_B)

  `CONSTANTES
C_REEL($vOflo_F;$vMinExp_F;$vMaxPrec_F)
$vMaxPrec_F:=16  ` Précision
$vOflo_F:=10^37
$vMinExp_F:=-87
$vIntegrale_F:=1
Si (($y_F>0) & ($ddl_F>0))
   `   Analyse la valeur de $vF_F
 $vF_F:=Log($y_F)
 $vF_F:=$ddl_F*$vF_F
 $vF_F:=$vF_F-MATH_LogGamma ($ddl_F+1)
 $vF_F:=$vF_F-$y_F
   ` soit $vF_F:=($ddl_F*Log($y_F))-MATH_LogGamma ($ddl_F+1)-$y_F
 Si ($vF_F>=$vMinExp_F)
  $vF_F:=Exp($vF_F)
  Si ($vF_F#0)
     ` Mise en place de la précision
   Au cas ou 
    : ($vDprec_I>$vMaxPrec_F)
    $vDprec_I:=$vMaxPrec_F
    : ($vDprec_I<=0)
    $vDprec_I:=1
   Fin de cas 
   $vCprec_F:=$vDprec_I
   $vEps_F:=10^-$vDprec_I
   
     `  Choisi entre les séries infinies ou les fractions continues
   Si (($y_F>1) & ($y_F>=$ddl_F))
      ` Fractions continues        
    $vA_F:=1-$ddl_F
    $B_F:=$vA_F+$y_F+1
    $vTerm_F:=0
    $tPn_F{1}:=1
    $tPn_F{2}:=$y_F
    $tPn_F{3}:=$y_F+1
    $tPn_F{4}:=$y_F*$B_F
    $vGin_F:=$tPn_F{3}/$tPn_F{4}
    $vDone_B:=Faux
    $vIter_I:=0
    Repeter 
     $vIter_I:=$vIter_I+1
     $vA_F:=$vA_F+1
     $B_F:=$B_F+2
     $vTerm_F:=$vTerm_F+1
     $vAn_F:=$vA_F*$vTerm_F
     $tPn_F{5}:=$B_F*$tPn_F{3}
     $tPn_F{5}:=$tPn_F{5}-($vAn_F*$tPn_F{1})
       ` soit $tPn_F{5}:=($B_F*$tPn_F{3})-($vAn_F*$tPn_F{1})
     $tPn_F{6}:=$B_F*$tPn_F{4}
     $tPn_F{6}:=$tPn_F{6}-($vAn_F*$tPn_F{2})
       ` Soit $tPn_F{6}:=($B_F*$tPn_F{4})-($vAn_F*$tPn_F{2})
     Si ($tPn_F{6}#0)
      $vRn_F:=$tPn_F{5}/$tPn_F{6}
      $vDif_F:=Abs($vGin_F-$vRn_F)
      Si ($vDif_F<=$vEps_F) & ($vDif_F<=($vEps_F*$vRn_F))
       $vDone_B:=Vrai
      Fin de si 
      $vGin_F:=$vRn_F
     Fin de si 
     $tPn_F{1}:=$tPn_F{3}
     $tPn_F{2}:=$tPn_F{4}
     $tPn_F{3}:=$tPn_F{5}
     $tPn_F{4}:=$tPn_F{6}
     Si (Abs($tPn_F{5})>=$vOflo_F)
      $tPn_F{1}:=$tPn_F{1}/$vOflo_F
      $tPn_F{2}:=$tPn_F{2}/$vOflo_F
      $tPn_F{3}:=$tPn_F{3}/$vOflo_F
      $tPn_F{4}:=$tPn_F{4}/$vOflo_F
     Fin de si 
    Jusque (($vIter_I>$vMaxIter_I) | $vDone_B)
    $vGin_F:=1-($vF_F*$vGin_F*$ddl_F)
    $vIntegrale_F:=$vGin_F
      `  Calcule la précision du résultat (débogage, non renvoyé)
      `
    Si ($vDif_F#0)
     $vCprec_F:=-(Log($vDif_F)/Log(10))  ` Ou MATH_LogN($vDif_F;10)
    Sinon 
     $vCprec_F:=$vMaxPrec_F
    Fin de si 
    
   Sinon 
      ` Séries infinies
    $vIter_I:=0
    $vTerm_F:=1
    $vC_F:=1
    $vA_F:=$ddl_F
    $vDone_B:=Faux
    Repeter 
     $vA_F:=$vA_F+1
     $vTerm_F:=($vTerm_F*$y_F)/$vA_F
     $vC_F:=$vC_F+$vTerm_F
     $vIter_I:=$vIter_I+1
    Jusque (($vIter_I>$vMaxIter_I) | (($vTerm_F/$vC_F)<=$vEps_F))
    $vIntegrale_F:=$vC_F*$vF_F
      `  Calcule la précision du résultat (débogage, non renvoyé)
    $vCprec_F:=-(Log($vTerm_F/$vC_F)/Log(10))  ` Ou MATH_LogN($vTerm_F/$vC_F;10)
   Fin de si 
   MATHERROR:=0  ` Tout s'est bien passé
  Sinon 
   MATHERROR:=-4  ` La valeur du paramètre intermédiaire exp($vF_F) est nulle
  Fin de si 
 Sinon 
  MATHERROR:=-3  ` La valeur du paramètre intermédiaire $vF_F est inférieure au minimum
 Fin de si 
Sinon 
 Si ($y_F<=0)
  MATHERROR:=-1  ` Valeur de la distribution négatif ou nul
 Sinon 
  MATHERROR:=-2  ` Degrés de liberté négatifs ou nuls
 Fin de si 
Fin de si 
$0:=$vIntegrale_F