From: N. Brouard Date: Wed, 15 Mar 2006 17:42:26 +0000 (+0000) Subject: (Module): Bug if status = -2, the loglikelihood was X-Git-Tag: imach-099s7~503 X-Git-Url: https://henry.ined.fr/git/?a=commitdiff_plain;h=79951640ff05ddd387dec193e57af25e084d70b7;p=.git (Module): Bug if status = -2, the loglikelihood was computed as likelihood omitting the logarithm. Version O.98e --- diff --git a/src/imach.c b/src/imach.c index f6dd1f9..13dcbac 100644 --- a/src/imach.c +++ b/src/imach.c @@ -1,6 +1,13 @@ /* $Id$ $State$ $Log$ + Revision 1.118 2006/03/14 18:20:07 brouard + (Module): varevsij Comments added explaining the second + table of variances if popbased=1 . + (Module): Covariances of eij, ekl added, graphs fixed, new html link. + (Module): Function pstamp added + (Module): Version 0.98d + Revision 1.117 2006/03/14 17:16:22 brouard (Module): varevsij Comments added explaining the second table of variances if popbased=1 . @@ -318,7 +325,7 @@ extern int errno; /* $Id$ */ /* $State$ */ -char version[]="Imach version 0.98d, March 2006, INED-EUROREVES-Institut de longevite "; +char version[]="Imach version 0.98e, March 2006, INED-EUROREVES-Institut de longevite "; char fullversion[]="$Revision$ $Date$"; char strstart[80]; char optionfilext[10], optionfilefiname[FILENAMELENGTH]; @@ -1376,20 +1383,20 @@ double func( double *x) } else if (s2==-2) { for (j=1,survp=0. ; j<=nlstate; j++) survp += out[s1][j]; - lli= survp; + lli= log(survp); } - else if (s2==-4) { - for (j=3,survp=0. ; j<=nlstate; j++) - survp += out[s1][j]; - lli= survp; - } +/* else if (s2==-4) { */ +/* for (j=3,survp=0. ; j<=nlstate; j++) */ +/* survp += out[s1][j]; */ +/* lli= survp; */ +/* } */ - else if (s2==-5) { - for (j=1,survp=0. ; j<=2; j++) - survp += out[s1][j]; - lli= survp; - } +/* else if (s2==-5) { */ +/* for (j=1,survp=0. ; j<=2; j++) */ +/* survp += out[s1][j]; */ +/* lli= survp; */ +/* } */ else{ @@ -1586,7 +1593,11 @@ double funcone( double *x) */ if( s2 > nlstate && (mle <5) ){ /* Jackson */ lli=log(out[s1][s2] - savm[s1][s2]); - } else if (mle==1){ + } else if (s2==-2) { + for (j=1,survp=0. ; j<=nlstate; j++) + survp += out[s1][j]; + lli= log(survp); + }else if (mle==1){ lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */ } else if(mle==2){ lli= (savm[s1][s2]>(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]):log((1.+bbh)*out[s1][s2])); /* linear interpolation */ @@ -1602,8 +1613,8 @@ double funcone( double *x) ll[s[mw[mi][i]][i]] += 2*weight[i]*lli; /* printf("i=%6d s1=%1d s2=%1d mi=%1d mw=%1d dh=%3d prob=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],out[s1][s2],savm[s1][s2]); */ if(globpr){ - fprintf(ficresilk,"%9d %6d %1d %1d %1d %1d %3d %10.6f %6.4f\ - %10.6f %10.6f %10.6f ", \ + fprintf(ficresilk,"%9d %6d %2d %2d %1d %1d %3d %11.6f %8.4f\ + %11.6f %11.6f %11.6f ", \ num[i],i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i], 2*weight[i]*lli,out[s1][s2],savm[s1][s2]); for(k=1,llt=0.,l=0.; k<=nlstate; k++){ @@ -2468,23 +2479,21 @@ void evsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, in for (age=bage; age<=fage; age ++){ /* If stepm=6 months */ /* Computed by stepm unit matrices, product of hstepm matrices, stored in an array of nhstepm length: nhstepm=10, hstepm=4, stepm=6 months */ - + hpxij(p3mat,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, cij); - + hf=hstepm*stepm/YEARM; /* Duration of hstepm expressed in year unit. */ - - /* Computing Variances of health expectancies */ - /* Gradient is computed with plus gp and minus gm. Code is duplicated in order to - decrease memory allocation */ - printf("%d|",(int)age);fflush(stdout); - fprintf(ficlog,"%d|",(int)age);fflush(ficlog); + + printf("%d|",(int)age);fflush(stdout); + fprintf(ficlog,"%d|",(int)age);fflush(ficlog); + /* Computing expectancies */ for(i=1; i<=nlstate;i++) for(j=1; j<=nlstate;j++) for (h=0, eij[i][j][(int)age]=0; h<=nhstepm-1; h++){ eij[i][j][(int)age] += (p3mat[i][j][h]+p3mat[i][j][h+1])/2.0*hf; -/* if((int)age==70)printf("i=%2d,j=%2d,h=%2d,age=%3d,%9.4f,%9.4f,%9.4f\n",i,j,h,(int)age,p3mat[i][j][h],hf,eij[i][j][(int)age]);*/ + /* if((int)age==70)printf("i=%2d,j=%2d,h=%2d,age=%3d,%9.4f,%9.4f,%9.4f\n",i,j,h,(int)age,p3mat[i][j][h],hf,eij[i][j][(int)age]);*/ } @@ -2498,12 +2507,12 @@ void evsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, in fprintf(ficreseij,"%9.4f", eip ); } fprintf(ficreseij,"\n"); - + } free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); printf("\n"); fprintf(ficlog,"\n"); - + } void cvevsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int cij, int estepm,double delti[],double **matcov,char strstart[] ) @@ -4404,7 +4413,7 @@ int main(int argc, char *argv[]) while ((val = strsep(&tok, "\"" )) != NULL && *val == '\0'); printf("val= |%s| pathr=%s\n",val,pathr); strcpy (pathtot, val); - if(pathr[0] == '\0') break; /* Un peu sale */ + if(pathr[0] == '\0') break; /* Dirty */ } } else{ @@ -4423,7 +4432,9 @@ int main(int argc, char *argv[]) /* Split argv[1]=pathtot, parameter file name to get path, optionfile, extension and name */ split(pathtot,path,optionfile,optionfilext,optionfilefiname); printf("\npathtot=%s,\npath=%s,\noptionfile=%s \noptionfilext=%s \noptionfilefiname=%s\n",pathtot,path,optionfile,optionfilext,optionfilefiname); - chdir(path); + chdir(path); /* Can be a relative path */ + if(getcwd(pathcd,MAXLINE) > 0) /* So pathcd is the full path */ + printf("Current directory %s!\n",pathcd); strcpy(command,"mkdir "); strcat(command,optionfilefiname); if((outcmd=system(command)) != 0){ @@ -5256,7 +5267,7 @@ Interval (in months) between two waves: Min=%d Max=%d Mean=%.2lf
\n",\ printf("%d %.0lf %lf %.0lf %.0lf %.0lf %lf\n",k,lsurv[k],p[1]*exp(p[2]*(k-agegomp)),(p[1]*exp(p[2]*(k-agegomp)))*lsurv[k],lpop[k],tpop[k],tpop[k]/lsurv[k]); - replace_back_to_slash(pathc,path); /* Even gnuplot wants a / */ + replace_back_to_slash(pathc,pathcd); /* Even gnuplot wants a / */ printinggnuplotmort(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p); printinghtmlmort(fileres,title,datafile, firstpass, lastpass, \ @@ -5493,7 +5504,7 @@ Interval (in months) between two waves: Min=%d Max=%d Mean=%.2lf
\n",\ /* freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint);*/ /*,dateprev1,dateprev2,jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/ - replace_back_to_slash(pathc,path); /* Even gnuplot wants a / */ + replace_back_to_slash(pathc,pathcd); /* Even gnuplot wants a / */ printinggnuplot(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p); printinghtml(fileres,title,datafile, firstpass, lastpass, stepm, weightopt,\ @@ -5888,7 +5899,12 @@ Interval (in months) between two waves: Min=%d Max=%d Mean=%.2lf
\n",\ fclose(ficlog); /*------ End -----------*/ - chdir(path); + + printf("Before Current directory %s!\n",pathcd); + if(chdir(pathcd) != 0) + printf("Can't move to directory %s!\n",path); + if(getcwd(pathcd,MAXLINE) > 0) + printf("Current directory %s!\n",pathcd); /*strcat(plotcmd,CHARSEPARATOR);*/ sprintf(plotcmd,"gnuplot"); #ifndef UNIX