suivant: Unité du modèle désagrégé
monter: Code source
précédent: Module du modèle gravitaire
Table des matières
// cette unité contient la définition de la strucuture nécessire
// au modèle et les procédures de manipulation des données dynamiques
// initialisées dans l'unité 'repartition'
unit repartition_utils;
interface
uses
Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms,
CheckLst, StdCtrls, Dialogs, Menus, Grids;
//------------------------------------------------------------------------------
//sructure de donnee dynamique associee aux parametres
//------------------------------------------------------------------------------
Type
cellule = RECORD
place : integer; //n° du champ dans la base [1..n]
utilite : integer; //1=UVP 2=UTC 3=UML
valeur : double; //valeur du parametre à estimer
modal : bool; // 1=variable a valeurs discretes 0 sinon
modal_val : integer; //valeur de la modalite correspondante
test : double;
end;
taille_modalite = integer;
Const
MAX = (65520 div Sizeof(cellule));
MAX_modalite = (65520 div Sizeof(taille_modalite));
Type
geniteur_model = ARRAY [1..MAX] of cellule;
geniteur_taille_modalite = ARRAY [1..MAX_modalite] of taille_modalite;
//------------------------------------------------------------------------------
//formulaire
//------------------------------------------------------------------------------
type
TForm4 = class(TForm)
MainMenu1: TMainMenu;
Fichier1: TMenuItem;
Sauver1: TMenuItem;
N1: TMenuItem;
Fermer1: TMenuItem;
StringGrid1: TStringGrid;
SaveDialog1: TSaveDialog;
procedure Fermer1Click(Sender: TObject);
procedure affiche(p : pointer);
procedure Sauver1Click(Sender: TObject);
private
{ Private declarations }
public
{ Public declarations }
end;
//------------------------------------------------------------------------------
//variables
//------------------------------------------------------------------------------
var
_place_mode : integer = 1;
nteta : integer;
nteta_vp : integer;
nteta_tc : integer;
nteta_ml : integer;
teta_view : TForm4;
teta : ^geniteur_model;
gk : ^geniteur_model;
mk : ^geniteur_model;
teta_temp : ^geniteur_model;
modalites : ^geniteur_taille_modalite;
lteta : double;
lzero : double;
compteur : integer;
//outils algebriques
procedure copie(source, cible : pointer);
procedure plus(cible, teta_1, teta_2 : pointer);
procedure fois(lambda : double; cible : pointer);
function norme(p : pointer) : double;
//outils fonctionnels
function pnj(p : pointer; u : integer) : double;
function f(p : pointer) : double;
function utilite(p : pointer; u : integer) : double;
procedure gradient(p : pointer);
implementation
uses repartition;
{$R *.DFM}
{******************************************************************************
* *
* *
* outils algebriques pour les tetas *
* *
* *
******************************************************************************}
procedure copie(source,cible : pointer);
var
i : integer;
s,c : ^geniteur_model;
begin
s:=source;
c:=cible;
reallocmem(s,nteta*sizeof(cellule));
reallocmem(c,nteta*sizeof(cellule));
for i := 1 to nteta do
begin
with c^[i] do
begin
valeur:=s^[i].valeur;
modal:=s^[i].modal;
modal_val:=s^[i].modal_val;
place:=s^[i].place;
utilite:=s^[i].utilite;
test:=s^[i].test;
end;
end;
end;
//-----------------------------------------------------------------------------
procedure plus(cible,teta_1,teta_2 : pointer);
var
i : integer;
c, t1, t2 : ^geniteur_model;
begin
c:=cible;
t1:=teta_1;
t2:=teta_2;
reallocmem(c,nteta*sizeof(cellule));
reallocmem(t1,nteta*sizeof(cellule));
reallocmem(t2,nteta*sizeof(cellule));
for i := 1 to nteta do c^[i].valeur:=t1^[i].valeur+t2^[i].valeur;
end;
//-----------------------------------------------------------------------------
procedure fois(lambda:double; cible : pointer);
var
i : integer;
c : ^geniteur_model;
begin
c:=cible;
reallocmem(c,nteta*sizeof(cellule));
for i := 1 to nteta do c^[i].valeur:=c^[i].valeur*lambda;
end;
//-----------------------------------------------------------------------------
// retourne la norme euclidienne d'un teta
//-----------------------------------------------------------------------------
function norme(p : pointer) : double;
var
i : integer;
temp : double;
x : ^geniteur_model;
begin
temp:=0;
x:=p;
reallocmem(x,nteta*sizeof(cellule));
for i := 1 to nteta do
begin
temp:=temp+sqr(x^[i].valeur);
end;
result:=sqrt(temp);
end;
{******************************************************************************
* *
* *
* outils d'analyse fonctionnelle pour l'algo du gradient conjugué *
* *
* *
******************************************************************************}
//-----------------------------------------------------------------------------
//calcule l'utilite pour un mode donné
//-----------------------------------------------------------------------------
function utilite(p : pointer; u : integer) : double;
var
i : integer;
temp : double;
x : ^geniteur_model;
begin
x:=p;
reallocmem(x,nteta*sizeof(cellule));
temp:=0;
//pour le mode VP
if (u=1) then
begin
temp:=temp+x^[1].valeur;
for i := 2 to nteta_vp do
begin
if x^[i].modal then
begin
if (x^[i].modal_val=Form3.echantillon.Fields[x^[i].place-1].Value) then
begin
temp:=temp+x^[i].valeur;
end;
end
else
begin
temp:=temp+x^[i].valeur*Form3.echantillon.Fields[x^[i].place-1].Value;
end;
end;
end;
//pour le mode TC
if (u=2) then
begin
temp:=temp+x^[1+nteta_vp].valeur;
for i := nteta_vp+2 to nteta_vp+nteta_tc do
begin
if x^[i].modal then
begin
if (x^[i].modal_val=Form3.echantillon.Fields[x^[i].place-1].Value) then
begin
temp:=temp+x^[i].valeur;
end;
end
else
begin
temp:=temp+x^[i].valeur*Form3.echantillon.Fields[x^[i].place-1].Value;
end;
end;
end;
//pour le ML
if (u=3) then
begin
for i := (nteta_tc+nteta_vp+1) to nteta do
begin
if x^[i].modal then
begin
if (x^[i].modal_val=Form3.echantillon.Fields[x^[i].place-1].Value) then
begin
temp:=temp+teta^[i].valeur;
end;
end
else
begin
temp:=temp+x^[i].valeur*Form3.echantillon.Fields[x^[i].place-1].Value;
end;
end;
end;
utilite:=temp;
end;
//-----------------------------------------------------------------------------
//calcule Pnj(teta,utilite)
//-----------------------------------------------------------------------------
function pnj(p : pointer; u : integer) : double;
var
temp : double;
uvp, utc, uml : double;
begin
temp:=0;
uvp:=utilite(p, 1);
utc:=utilite(p, 2);
uml:=utilite(p, 3);
case u of
1 : temp := uvp;
2 : temp := utc;
3 : temp := uml;
end;
// exp(u_mode_choisi)
// renvoie ------------------
// somme exp(u_modes)
pnj:=exp(temp)/(exp(uvp)+exp(utc)+exp(uml));
end;
//-----------------------------------------------------------------------------
//calcule la vraisemblance de l'échantillon pour un teta donné.
//-----------------------------------------------------------------------------
function f(p : pointer) : double;
var
i : integer;
u : integer;
vraisemblance : double;
temp : double;
begin
Form3.echantillon.First;
vraisemblance:=0;
for i := 1 to Form3.echantillon.RecordCount do
begin
u:=Form3.echantillon.Fields[_place_mode-1].value;
temp:=ln(pnj(p,u));
vraisemblance:=vraisemblance+temp;
Form3.echantillon.Next;
end;
result:=vraisemblance;
end;
//-----------------------------------------------------------------------------
// calcule le gradient de la vraisemblance pour un teta donné
//-----------------------------------------------------------------------------
procedure gradient(p : pointer);
var
test : boolean;//=Vrai si i_teta est associe a l'utilite u de l'individu n
n : integer;//individu actif
i_teta : integer;//compteur pour les coordonnees de teta.
u : integer;//utilite coisie par n
p1,p2 : double;//les deux parties de la derivee partielle
x : ^geniteur_model;
xink : integer;
begin
x:=p;
reallocmem(x,nteta*sizeof(cellule));
xink:=0;
// il faut initialiser gk
for i_teta := 1 to nteta do gk^[i_teta].valeur:=0;
// calculer chaque composante du gradient ...
for i_teta := 1 to nteta do
begin
Form3.echantillon.First;
// ... en parcourant la base
for n := 1 to Form3.echantillon.RecordCount do
begin
u:=Form3.echantillon.Fields[_place_mode-1].value;
//activer la variable test si le teta actif est ds l'utilite du mode
// que n a choisi
if u=x^[i_teta].utilite then test:=True else test:=False;
//xink
if x^[i_teta].modal then
begin
if x^[i_teta].modal_val=Form3.echantillon.Fields[x^[i_teta].place-1].Value
then xink:=1
else xink:=0;
end;
//calcul de p1
if test then
begin
if ((i_teta=1) or (i_teta=nteta_vp+1)) then p1:=1
else
if x^[i_teta].modal then p1:=xink
else p1:=Form3.echantillon.Fields[x^[i_teta].place-1].Value;
end
else p1:=0;
//calcul de p2
if ((i_teta=1) or (i_teta=nteta_vp+1))
then p2:=exp(utilite(x,x^[i_teta].utilite))
else
begin
if x^[i_teta].modal
then p2:=xink*exp(utilite(x,x^[i_teta].utilite))
else p2:=Form3.echantillon.Fields[x^[i_teta].place-1].Value
*exp(utilite(x,x^[i_teta].utilite));
end;
p2:=p2/(exp(utilite(x,1))+exp(utilite(x,2))+exp(utilite(x,3)));
gk^[i_teta].valeur:=gk^[i_teta].valeur+p1-p2;
Form3.echantillon.Next;
end;
end;
end;
//*****************************************************************************
//*****************************************************************************
//*****************************************************************************
procedure TForm4.Fermer1Click(Sender: TObject);
begin
Close;
end;
//-----------------------------------------------------------------------------
procedure TForm4.affiche(p : pointer);
var
i : integer;
x : ^geniteur_model;
begin
x:=p;
reallocmem(x,nteta*sizeof(cellule));
StringGrid1.RowCount:=nteta+3;
StringGrid1.Colcount:=7;
StringGrid1.Cells[1,0]:='valeur';
StringGrid1.Cells[2,0]:='n° du Champ';
StringGrid1.Cells[3,0]:='utilité';
StringGrid1.Cells[4,0]:='variable modale';
StringGrid1.Cells[5,0]:='valeur de la modalité';
StringGrid1.Cells[6,0]:='test';
for i := 1 to nteta do
begin
StringGrid1.Cells[0,i]:='teta '+inttostr(i);
StringGrid1.Cells[1,i]:=floattostr(x^[i].valeur);
if x^[i].place=0
then StringGrid1.Cells[2,i]:='constante'
else
StringGrid1.Cells[2,i]:=Form3.echantillon.Fields[x^[i].place-1].FullName;
case x^[i].utilite of
1 : StringGrid1.Cells[3,i]:='vp';
2 : StringGrid1.Cells[3,i]:='tc';
3 : StringGrid1.Cells[3,i]:='ml';
end;
if x^[i].modal then
begin
StringGrid1.Cells[4,i]:='oui';
StringGrid1.Cells[5,i]:=inttostr(x^[i].modal_val);
end
else StringGrid1.Cells[4,i]:='non';
StringGrid1.Cells[6,i]:=floattostr(x^[i].test);
end;
StringGrid1.Cells[1,nteta+1]:='vraisemblance';
StringGrid1.Cells[2,nteta+1]:='vraisemblance(0)';
StringGrid1.Cells[3,nteta+1]:='rho2';
StringGrid1.Cells[4,nteta+1]:='itérations';
StringGrid1.Cells[1,nteta+2]:=floattostr(lteta);
StringGrid1.Cells[2,nteta+2]:=floattostr(lzero);
StringGrid1.Cells[3,nteta+2]:=floattostr(1-lteta/lzero);
StringGrid1.Cells[4,nteta+2]:=floattostr(compteur);
end;
//------------------------------------------------------------------------------
procedure TForm4.Sauver1Click(Sender: TObject);
var
fichier : textfile;
i,j : integer;
begin
if savedialog1.execute then
begin
AssignFile(fichier,savedialog1.FileName);
Rewrite(fichier);
for i := 1 to Stringgrid1.ColCount+1 do
begin
for j := 1 to Stringgrid1.RowCount do write(fichier,stringgrid1.cells[j-1,i-1]+' ');
writeln(fichier,'');
end;
closefile(fichier);
end;
end;
//------------------------------------------------------------------------------
end.
2003-06-21