Vous n'êtes pas identifié.
Vous êtes sur le forum du master ESA !
Le site du master ESA - description de la formation, notes de cours, contacts... vient de déménager !!!
Venez visiter notre nouveau site : www.master-esa.fr
D'abord les conditions initiales issues d'un logit :
%macro paramin(base=,depvar=,listvar=); proc logistic data=&base; ods output ParameterEstimates=parametres; model &depvar = &listvar ; run; data param (keep = estimate);set parametres;run; proc transpose data=param out=paramin;run; data paramin (drop= _NAME_ _LABEL_); set paramin;run; %mend; %paramin(base=db.insample,depvar=incident,listvar=revenu duree montcred ancienn);
ensuite l'estimation semi-parametrique :
%macro klein_spady(base=,base2=,kernel=,listvar=,h=); proc iml; use &base; read all var { &listvar } into x; use &base2; read all var _all_ into parmin; start loglik(theta) global(x); nbrlin=nrow(x); nbrcol=ncol(x); ll=0; m=j(nbrlin,3,0); do i=1 to nbrlin; bidul=0; do k=1 to nbrcol-1; bidulk=theta[k+1]*x[i,k]; bidul=bidul+bidulk; end; scorei=bidul+theta[1]; wn=0; wd=0; do j=1 to nbrlin; machin=0; do l=1 to nbrcol-1; machinl=theta[l+1]*x[j,l]; machin=machin+machinl; end; scorej=machin+theta[1]; kervar=(scorej-scorei)/&h; %if &kernel=normal %then %do; kernelj=(1/sqrt(2*3.14))*exp(-0.5*(kervar)##2); %end; %if &kernel=epanechnicov %then %do; if -1<=kervar<=1 then do; kernelj=3/4*(1-kervar##2); end; else do; kernelj=0; end; %end; %if &kernel=triangulaire %then %do; if -1<=kervar<=1 then do; kernelj=1-abs(kervar); end; else do; kernelj=0; end; %end; %if &kernel=uniforme %then %do; if -1<=kervar<=1 then do; kernelj=0.5; end; else do; kernelj=0; end; %end; w1=kernelj*x[j,nbrcol]; wn=wn+w1; wd=wd+kernelj; end; ri=wn/wd; m[i,1]=ri; m[i,2]=scorei; m[i,3]=x[i,nbrcol]; lli=x[i,nbrcol]*log(ri)+(1-x[i,nbrcol])*log(1-ri); ll=ll+lli; end; create sortie FROM m ; append FROM m ; close sortie ; return(ll); finish loglik; options = {1,2}; call nlpnra(rc,thetahat,"loglik",parmin,options); print," Résultats de l estimation par Maximum de Vraisemblance "; quit; %mend; %klein_spady(base=db.insample,base2=paramin,kernel=epanechnicov,listvar=revenu duree montcred ancienn age incident,h=0.35);
Pour tracer la fonction de répartition estimée:
proc sort data=sortie; by score; run; symbol i=none value=diamond height=0.75; axis1 label=('Score' justify=right); axis2 label=('Proba' justify=right); proc gplot data=sortie; plot COL1*COL2 /vref=0.5 haxis=axis1 vaxis=axis2; title 'Fonction de répartition'; run; quit;
CAUTION : Ce n'est pas encore corrigé ! donc je ne sais pas si tout est bon dans la macro ! je la poste pour ceux parmis vous qui luttent sur leurs memoires lol ( j'ai pas vos mail ). Pour la roc curve il y a les probas estimées dans la table 'sortie' COL1 ( soyez cool j'ai du taff !!!).
Hors ligne
Petite correction :
l'instruction dans le module « loglik »
create sortie FROM m ;
append FROM m ;
close sortie ;
NE SERT A RIEN ! Désolé !
Faites plutôt un " RUN loglik(vos_parametres) ; " puis recréer la table sortie avec l'étape précédente, une fois que vous avez vos paramètres estimés bien sûr...
Ça marche comment Paypal ???
Hors ligne
Hors ligne
LOL !
on apprend beaucoup de chose sur vous en regardant d'où viennent les images que vous intégrez à vos messages :
http://www.vaq.qc.ca/
voitures anciennes du quebec !
cordialement
SR
Hors ligne
Et oui ! fan du PSG... Collectionneur de vielles bagnoles québécoises...ça fait pas mal de defaut tout ca!!
Je vais commencer à verifier les liens avant de poster mes images lol
Hors ligne
alaa-eddine a écrit:
Et oui ! fan du PSG...
grrrrrrrrrrrrrrrr beurk
est ce une equipe ca!!!! le psg et oui paris n'est plus magique je dirai plutot paris c'est tragique!!!! vive l 'om
bon je me remets au boulot lol...
Hors ligne
MAIS NON !!! Je ne suis pas fan du PSG !!!!! ni de l'OM d'ailleurs et encore moins de Lyon...lol et bien sûr loin d'être Collectionneur de voitures québécoises ! je ne suis pas encore à la retraite...
Hors ligne
success story ! lol
Hors ligne
alaa-eddine a écrit:
Petite correction :
l'instruction dans le module « loglik »
create sortie FROM m ;
append FROM m ;
close sortie ;
NE SERT A RIEN ! Désolé !
Faites plutôt un " RUN loglik(vos_parametres) ; " puis recréer la table sortie avec l'étape précédente, une fois que vous avez vos paramètres estimés bien sûr...
Ça marche comment Paypal ???
Salut,
Je viens de voir votre programme. Ce que vous écrivez ci-dessus me parait confus. Vous pouvez être un peu plus explicite svp?
Merci
Hors ligne
En fait le but de l'instruction :
create sortie FROM m ;
append FROM m ;
close sortie ;
est de récupérer les r_chapeau, les x_beta et la variable expliquée.
Et donc pour aller plus vite (lorsque la taille de l'échantillon est importante), je propose dans un premier temps de faire tourner la macro sans cette instruction, pour estimer les paramètres.
Ensuite une fois que vous avez vos paramètres estimés vous refaite tourner la macro en mettant comme argument au module "loglik" ces paramètres estimés.
Voici les dernières ligne de la macro (imaginons que mes paramètres estimés sont : -0.347016 0.408833 -0.336351 -0.431585 0.015167 -0.346089 -0.014785 , donc 7 régresseurs dont la constante ):
...
ri=wn/wd;
m[i,1]=ri;
m[i,2]=scorei;
m[i,3]=x[i,nbrcol];
create sortie FROM m ;
append FROM m ;
close sortie ;
finish loglik;
parametres_est={-0.347016 0.408833 -0.336351 -0.431585 0.015167 -0.346089 -0.014785};
RUN loglik(parametres_est) ; *au lieu du call nlpnra(rc,thetahat,"loglik",parmin,options);
quit;
%mend;
J'espère que c'est plus clair maintenant...
Hors ligne
Il manque juste un "end;" :
...
ri=wn/wd;
m[i,1]=ri;
m[i,2]=scorei;
m[i,3]=x[i,nbrcol];
end;
create sortie FROM m ;
append FROM m ;
close sortie ;
finish loglik;
parametres_est={-0.347016 0.408833 -0.336351 -0.431585 0.015167 -0.346089 -0.014785};
RUN loglik(parametres_est) ; *au lieu du call nlpnra(rc,thetahat,"loglik",parmin,options);
quit;
%mend;
C'est good ?
Hors ligne
alaa-eddine a écrit:
parametres_est={-0.347016 0.408833 -0.336351 -0.431585 0.015167 -0.346089 -0.014785};
RUN loglik(parametres_est) ; *au lieu du call nlpnra(rc,thetahat,"loglik",parmin,options);
quit;
%mend;
C'est good ?
Cher Alaa-eddine
1. si tu as 30 régresseurs, tu les réecris tous dans "parametres-est="?
2. as tu une idée du délai moyen de l'estimation (30000 observations)?
Merci bcp
Mak
Hors ligne
Re tout le monde !!!
Salut Mak ! Désolé pour le retard !
Alors la solution est simple !
Il suffit de rajouter cette ligne de commande pour récupérer le vecteur des « beta chapeau » (après le « call nlpqn »), ce qui donne :
call nlpnra(rc, betahat,"loglik",parmin,options);
call nlpfdd(f,g,hess,"loglik", betahat);
“betahat” sera le vecteur des 30 beta chapeau…
Ensuite tu les exporte de la proc iml :
create beta FROM betahat;
append FROM betahat;
close beta ;
il suffit maintenant de réimporter ces estimateur dans la 2ème proc iml :
use beta ;
read all var _all_ into betac ;
…
Et hop ! :
RUN loglik(betac) ;
quit;
%mend;
Pour les 30000 c’est tendu ! pour 2 raisons :
Faut savoir que le module iml est super instable lorsqu’il s’agit de faire trop de calcul, surtout si t’as 30 régresseurs… lol.
Tu risque de rejeter H0 tout le temps.
Je ferais plutôt un tirage aléatoire de 1500 ou 2000 individus maxi, pour ton échantillon insample, puis un échantillon test et/ou validation.
Pas facile à utiliser cette macro… surtout qu'il y a moyen de l'améliorer !
Si j’avais le temps je la referais en data step… mais bon… tendu en ce moment !
Hors ligne
Salut Alaa-Eddine,
Merci c bcp plus clair. En effet pas facile cette macro. Je la teste et te tiens au courant. Au fait c quoi la 2ème proc IML. C celle qui correspond à un kernel non gaussien?
il faut écrire thetahat au lieu de betahat. Dans le prog, tu as:
proc IML;
use &base;
read all var { &listvar } into x;
use &base2;
read all var _all_ into parmin;
start loglik(theta) global(x);
nbrlin=nrow(x);
nbrcol=ncol(x);
ll=0;
m=j(nbrlin,3,0);
do i=1 to nbrlin;
bidul=0;
do k=1 to nbrcol-1;
bidulk=theta[k+1]*x[i,k];
bidul=bidul+bidulk;
end;
scorei=bidul+theta[1];....
Dernière modification par mak (09-09-2008 01:57:51)
Hors ligne
J'ai essayé de faire tourner ton prog et voilà le message du "log" de sas après les dernières modifications que tu as proposées. qu'en penses tu?
NOTE: IML Ready
NOTE: I/O required temporary file to be opened.
NOTE: Module LOGLIK defined.
NOTE: The data set WORK.SORTIE has 3496 observations and 3 variables.
ERROR: Execution error as noted previously. (rc=100)
operation : NLPNRA at line 751 column 1
operands : *LIT1031, PARMIN, OPTIONS
*LIT1031 1 row 1 col (character, size 6)
loglik
PARMIN 1 row 49 cols (numeric)
OPTIONS 2 rows 1 col (numeric)
1
2
statement : CALL at line 751 column 1
ERROR: (execution) Matrix has not been set to a value.
operation : NLPFDD at line 751 column 1
operands : *LIT1032, THETAHAT
*LIT1032 1 row 1 col (character, size 6)
loglik
THETAHAT 0 row 0 col (type ?, size 0)
statement : CALL at line 751 column 1
ERROR: Matrix THETAHAT has not been set to a value.
statement : CREATE at line 751 column 1
ERROR: No data set is currently open for output.
statement : APPEND at line 751 column 1
NOTE: Cannot close WORK.BETA; it is not open.
Hors ligne
Bon, j'arrête de te saouler mais pourrais-tu une fois pour toute réecrire ton prog? je l'ai modifié pas mal de fois au point où je ne c plus ce qui ne va pas. Merci bcp m'sieur.
Hors ligne
ahhhh si j'avais le temps...
je ne te promet rien, car j'ai des semaines ultra chargées, mais je vais essayer !
juste une question : pourquoi du semi-paramétrique ? ça se trouve, un simple logit fait mieux... surtout si la macro ne fonctionne pas de manière optimale
Hors ligne
Double soucis des modèles dynamiques à la Heckman: 1. Conditions initiales et par csq modèle en équilibre douteux. 2. présence de corrélations. En panel, si T est grand (cf Hsiao, 1986, 2003) ces 2 points disparaissent (presque). G un T petit d'où solution modèle semi-paramétrique (cf Vella, Wooldridge). J'espère avoir répondu à ta question et merci pour on aide.
Hors ligne
SOS MACRO ALLAHEDHINE,ON COMPTE SUR TOI...la macro marche pas mak a raison ,ça me file le meme sms..
Dernière modification par boualdjasonia (13-04-2009 23:19:32)
Hors ligne
J’ai plus SAS à la maison (licence sas academic expirée)... faut que j'ai le temps de me replonger dedans au bureau ce qui est loin d’être gagné, i'll do my best !
Hors ligne