--- imach/src/imach.c 2015/10/23 15:50:53 1.205 +++ imach/src/imach.c 2015/10/27 17:36:57 1.207 @@ -1,6 +1,12 @@ -/* $Id: imach.c,v 1.205 2015/10/23 15:50:53 brouard Exp $ +/* $Id: imach.c,v 1.207 2015/10/27 17:36:57 brouard Exp $ $State: Exp $ $Log: imach.c,v $ + Revision 1.207 2015/10/27 17:36:57 brouard + *** empty log message *** + + Revision 1.206 2015/10/24 07:14:11 brouard + *** empty log message *** + Revision 1.205 2015/10/23 15:50:53 brouard Summary: 0.98r3 some clarification for graphs on likelihood contributions @@ -751,12 +757,12 @@ typedef struct { #define ODIRSEPARATOR '\\' #endif -/* $Id: imach.c,v 1.205 2015/10/23 15:50:53 brouard Exp $ */ +/* $Id: imach.c,v 1.207 2015/10/27 17:36:57 brouard Exp $ */ /* $State: Exp $ */ #include "version.h" char version[]=__IMACH_VERSION__; char copyright[]="October 2015,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015"; -char fullversion[]="$Revision: 1.205 $ $Date: 2015/10/23 15:50:53 $"; +char fullversion[]="$Revision: 1.207 $ $Date: 2015/10/27 17:36:57 $"; char strstart[80]; char optionfilext[10], optionfilefiname[FILENAMELENGTH]; int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */ @@ -1960,7 +1966,24 @@ double **prevalim(double **prlim, int nl { /* Computes the prevalence limit in each live state at age x by left multiplying the unit matrix by transitions matrix until convergence is reached with precision ftolpl */ - + /* Wx= Wx-1 Px-1= Wx-2 Px-2 Px-1 = Wx-n Px-n ... Px-2 Px-1 I */ + /* Wx is row vector: population in state 1, population in state 2, population dead */ + /* or prevalence in state 1, prevalence in state 2, 0 */ + /* newm is the matrix after multiplications, its rows are identical at a factor */ + /* Initial matrix pimij */ + /* {0.85204250825084937, 0.13044499163996345, 0.017512500109187184, */ + /* 0.090851990222114765, 0.88271245433047185, 0.026435555447413338, */ + /* 0, 0 , 1} */ + /* + * and after some iteration: */ + /* {0.45504275246439968, 0.42731458730878791, 0.11764266022681241, */ + /* 0.45201005341706885, 0.42865420071559901, 0.11933574586733192, */ + /* 0, 0 , 1} */ + /* And prevalence by suppressing the deaths are close to identical rows in prlim: */ + /* {0.51571254859325999, 0.4842874514067399, */ + /* 0.51326036147820708, 0.48673963852179264} */ + /* If we start from prlim again, prlim tends to a constant matrix */ + int i, ii,j,k; double min, max, maxmin, maxmax,sumnew=0.; /* double **matprod2(); */ /* test */ @@ -2016,13 +2039,13 @@ double **prevalim(double **prlim, int nl prlim[i][j]= newm[i][j]/(1-sumnew); max=FMAX(max,prlim[i][j]); min=FMIN(min,prlim[i][j]); - /* printf(" age= %d prevalim i=%d, j=%d, prmlim[%d][%d]=%f, agefin=%d max=%f min=%f\n", (int)age, i, j, i, j, prlim[i][j],(int)agefin, max, min); */ + printf(" age= %d prevalim i=%d, j=%d, prmlim[%d][%d]=%f, agefin=%d max=%f min=%f\n", (int)age, i, j, i, j, prlim[i][j],(int)agefin, max, min); } maxmin=(max-min)/(max+min)*2; maxmax=FMAX(maxmax,maxmin); } /* j loop */ *ncvyear= (int)age- (int)agefin; - /* printf("maxmax=%lf maxmin=%lf ncvloop=%ld, age=%d, agefin=%d ncvyear=%d \n", maxmax, maxmin, ncvloop, (int)age, (int)agefin, *ncvyear); */ + printf("maxmax=%lf maxmin=%lf ncvloop=%d, age=%d, agefin=%d ncvyear=%d \n", maxmax, maxmin, ncvloop, (int)age, (int)agefin, *ncvyear); if(maxmax < ftolpl){ /* printf("maxmax=%lf maxmin=%lf ncvloop=%ld, age=%d, agefin=%d ncvyear=%d \n", maxmax, maxmin, ncvloop, (int)age, (int)agefin, *ncvyear); */ return prlim; @@ -2647,7 +2670,7 @@ void likelione(FILE *ficres,double p[], fprintf(ficlog,"Problem with resultfile: %s\n", fileresilk); } fprintf(ficresilk, "#individual(line's_record) count age s1 s2 wave# effective_wave# number_of_matrices_product pij weight weight/gpw -2ln(pij)*weight 0pij_x 0pij_(x-stepm) cumulating_loglikeli_by_health_state(reweighted=-2ll*weightXnumber_of_contribs/sum_of_weights) and_total\n"); - fprintf(ficresilk, "#num_i age i s1 s2 mi mw dh likeli weight %weight 2wlli out sav "); + fprintf(ficresilk, "#num_i age i s1 s2 mi mw dh likeli weight %%weight 2wlli out sav "); /* 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; k<=nlstate; k++) fprintf(ficresilk," -2*gipw/gsw*weight*ll[%d]++",k); @@ -2662,17 +2685,17 @@ void likelione(FILE *ficres,double p[], else if(mle >=1) fprintf(fichtm,"\n
File of contributions to the likelihood computed with optimized parameters mle = %d.",mle); fprintf(fichtm," You should at least run with mle >= 1 to get starting values corresponding to the optimized parameters in order to visualize the real contribution of each individual/wave: %s
\n",subdirf(fileresilk),subdirf(fileresilk)); - - fprintf(fichtm,"
- The function drawn is -2Log(L) in Log scale: by state of origin %s-ori.png
\ + + fprintf(fichtm,"
- The function drawn is -2Log(L) in Log scale: by state of origin %s-ori.png
\ ",subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_")); - fprintf(fichtm,"
- and by state of destination %s-dest.png
\ + fprintf(fichtm,"
- and by state of destination %s-dest.png
\ ",subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_")); - fflush(fichtm); + fflush(fichtm); - for (k=1; k<= nlstate ; k++) { - fprintf(fichtm,"
- Probability p%dj by origin %d and destination j %s-p%dj.png
\ + for (k=1; k<= nlstate ; k++) { + fprintf(fichtm,"
- Probability p%dj by origin %d and destination j %s-p%dj.png
\ ",k,k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k); - } + } } return; }