--- imach/src/imach.c 2017/04/26 06:01:29 1.264 +++ imach/src/imach.c 2017/04/26 16:22:11 1.265 @@ -1,6 +1,9 @@ -/* $Id: imach.c,v 1.264 2017/04/26 06:01:29 brouard Exp $ +/* $Id: imach.c,v 1.265 2017/04/26 16:22:11 brouard Exp $ $State: Exp $ $Log: imach.c,v $ + Revision 1.265 2017/04/26 16:22:11 brouard + Summary: imach 0.99r13 Some bugs fixed + Revision 1.264 2017/04/26 06:01:29 brouard Summary: Labels in graphs @@ -994,12 +997,12 @@ typedef struct { #define ODIRSEPARATOR '\\' #endif -/* $Id: imach.c,v 1.264 2017/04/26 06:01:29 brouard Exp $ */ +/* $Id: imach.c,v 1.265 2017/04/26 16:22:11 brouard Exp $ */ /* $State: Exp $ */ #include "version.h" char version[]=__IMACH_VERSION__; char copyright[]="February 2016,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2018"; -char fullversion[]="$Revision: 1.264 $ $Date: 2017/04/26 06:01:29 $"; +char fullversion[]="$Revision: 1.265 $ $Date: 2017/04/26 16:22:11 $"; char strstart[80]; char optionfilext[10], optionfilefiname[FILENAMELENGTH]; int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */ @@ -4308,7 +4311,7 @@ void freqsummary(char fileres[], double int firstpass, int lastpass, int stepm, int weightopt, char model[]) { /* Some frequencies as well as proposing some starting values */ - int i, m, jk, j1, bool, z1,j, nj, nl, k, iv, jj=0; + int i, m, jk, j1, bool, z1,j, nj, nl, k, iv, jj=0, s1=1, s2=1; int iind=0, iage=0; int mi; /* Effective wave */ int first; @@ -4387,23 +4390,48 @@ Title=%s
Datafile=%s Firstpass=%d La k2cpt=0; if(cptcoveff == 0 ) - nl=1; /* Constant model only */ + nl=1; /* Constant and age model only */ else nl=2; + + /* if a constant only model, one pass to compute frequency tables and to write it on ficresp */ + /* Loop on nj=1 or 2 if dummy covariates j!=0 + * Loop on j1(1 to 2**cptcoveff) covariate combination + * freq[s1][s2][iage] =0. + * Loop on iind + * ++freq[s1][s2][iage] weighted + * end iind + * if covariate and j!0 + * headers Variable on one line + * endif cov j!=0 + * header of frequency table by age + * Loop on age + * pp[s1]+=freq[s1][s2][iage] weighted + * pos+=freq[s1][s2][iage] weighted + * Loop on s1 initial state + * fprintf(ficresp + * end s1 + * end age + * if j!=0 computes starting values + * end compute starting values + * end j1 + * end nl + */ for (nj = 1; nj <= nl; nj++){ /* nj= 1 constant model, nl number of loops. */ if(nj==1) j=0; /* First pass for the constant */ - else + else{ j=cptcoveff; /* Other passes for the covariate values */ + } first=1; - for (j1 = 1; j1 <= (int) pow(2,j); j1++){ /* Loop on covariates combination in order of model, excluding quantitatives, V4=0, V3=0 for example, fixed or varying covariates */ + for (j1 = 1; j1 <= (int) pow(2,j); j1++){ /* Loop on all covariates combination of the model, excluding quantitatives, V4=0, V3=0 for example, fixed or varying covariates */ posproptt=0.; /*printf("cptcoveff=%d Tvaraff=%d", cptcoveff,Tvaraff[1]); scanf("%d", i);*/ for (i=-5; i<=nlstate+ndeath; i++) - for (jk=-5; jk<=nlstate+ndeath; jk++) + for (s2=-5; s2<=nlstate+ndeath; s2++) for(m=iagemin; m <= iagemax+3; m++) - freq[i][jk][m]=0; + freq[i][s2][m]=0; for (i=1; i<=nlstate; i++) { for(m=iagemin; m <= iagemax+3; m++) @@ -4421,7 +4449,7 @@ Title=%s
Datafile=%s Firstpass=%d La /* dateintsum=0; */ /* k2cpt=0; */ - /* For that combination of covariate j1, we count and print the frequencies in one pass */ + /* For that combination of covariates j1 (V4=1 V3=0 for example), we count and print the frequencies in one pass */ for (iind=1; iind<=imx; iind++) { /* For each individual iind */ bool=1; if(j !=0){ @@ -4437,7 +4465,7 @@ Title=%s
Datafile=%s Firstpass=%d La /* /\* sumnew+=coqvar[z1][iind]; *\/ */ /* }else */ if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]){ /* for combination j1 of covariates */ - /* Tests if this individual iind responded to combination j1 (V4=1 V3=0) */ + /* Tests if the value of the covariate z1 for this individual iind responded to combination j1 (V4=1 V3=0) */ bool=0; /* bool should be equal to 1 to be selected, one covariate value failed */ /* printf("bool=%d i=%d, z1=%d, Tvaraff[%d]=%d, covar[Tvarff][%d]=%2f, codtabm(%d,%d)=%d, nbcode[Tvaraff][codtabm(%d,%d)=%d, j1=%d\n", bool,i,z1, z1, Tvaraff[z1],i,covar[Tvaraff[z1]][i],j1,z1,codtabm(j1,z1), @@ -4448,7 +4476,7 @@ Title=%s
Datafile=%s Firstpass=%d La } /* cptcovn > 0 */ } /* end any */ }/* end j==0 */ - if (bool==1){ /* We selected an individual iind satisfying combination j1 or all fixed */ + if (bool==1){ /* We selected an individual iind satisfying combination j1 (V4=1 V3=0) or all fixed covariates */ /* for(m=firstpass; m<=lastpass; m++){ */ for(mi=1; miDatafile=%s Firstpass=%d La /* fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/ - pstamp(ficresp); + if(cptcoveff==0 && nj==1) /* no covariate and first pass */ + pstamp(ficresp); if (cptcoveff>0 && j!=0){ + pstamp(ficresp); printf( "\n#********** Variable "); fprintf(ficresp, "\n#********** Variable "); fprintf(ficresphtm, "\n

********** Variable "); @@ -4539,20 +4569,23 @@ Title=%s
Datafile=%s Firstpass=%d La fprintf(ficlog, "**********\n"); } fprintf(ficresphtm,""); + if((cptcoveff==0 && nj==1)|| nj==2 ) /* no covariate and first pass */ + fprintf(ficresp, " Age"); + if(nj==2) for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresp, " V%d=%d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); for(i=1; i<=nlstate;i++) { - fprintf(ficresp, " Age Prev(%d) N(%d) N ",i,i); + if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp," Prev(%d) N(%d) N ",i,i); fprintf(ficresphtm, "",i,i); } - fprintf(ficresp, "\n"); + if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp, "\n"); fprintf(ficresphtm, "\n"); /* Header of frequency table by age */ fprintf(ficresphtmfr,"
AgePrev(%d)N(%d)N
"); fprintf(ficresphtmfr," "); - for(jk=-1; jk <=nlstate+ndeath; jk++){ + for(s2=-1; s2 <=nlstate+ndeath; s2++){ for(m=-1; m <=nlstate+ndeath; m++){ - if(jk!=0 && m!=0) - fprintf(ficresphtmfr," ",jk,m); + if(s2!=0 && m!=0) + fprintf(ficresphtmfr," ",s2,m); } } fprintf(ficresphtmfr, "\n"); @@ -4577,97 +4610,112 @@ Title=%s
Datafile=%s Firstpass=%d La fprintf(ficresphtmfr," ",iage); fprintf(ficlog,"Age %d", iage); } - for(jk=1; jk <=nlstate ; jk++){ - for(m=-1, pp[jk]=0; m <=nlstate+ndeath ; m++) - pp[jk] += freq[jk][m][iage]; + for(s1=1; s1 <=nlstate ; s1++){ + for(m=-1, pp[s1]=0; m <=nlstate+ndeath ; m++) + pp[s1] += freq[s1][m][iage]; } - for(jk=1; jk <=nlstate ; jk++){ + for(s1=1; s1 <=nlstate ; s1++){ for(m=-1, pos=0; m <=0 ; m++) - pos += freq[jk][m][iage]; - if(pp[jk]>=1.e-10){ + pos += freq[s1][m][iage]; + if(pp[s1]>=1.e-10){ if(first==1){ - printf(" %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]); + printf(" %d.=%.0f loss[%d]=%.1f%%",s1,pp[s1],s1,100*pos/pp[s1]); } - fprintf(ficlog," %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]); + fprintf(ficlog," %d.=%.0f loss[%d]=%.1f%%",s1,pp[s1],s1,100*pos/pp[s1]); }else{ if(first==1) - printf(" %d.=%.0f loss[%d]=NaNQ%%",jk,pp[jk],jk); - fprintf(ficlog," %d.=%.0f loss[%d]=NaNQ%%",jk,pp[jk],jk); + printf(" %d.=%.0f loss[%d]=NaNQ%%",s1,pp[s1],s1); + fprintf(ficlog," %d.=%.0f loss[%d]=NaNQ%%",s1,pp[s1],s1); } } - for(jk=1; jk <=nlstate ; jk++){ - /* posprop[jk]=0; */ - for(m=0, pp[jk]=0; m <=nlstate+ndeath; m++)/* Summing on all ages */ - pp[jk] += freq[jk][m][iage]; - } /* pp[jk] is the total number of transitions starting from state jk and any ending status until this age */ + for(s1=1; s1 <=nlstate ; s1++){ + /* posprop[s1]=0; */ + for(m=0, pp[s1]=0; m <=nlstate+ndeath; m++)/* Summing on all ages */ + pp[s1] += freq[s1][m][iage]; + } /* pp[s1] is the total number of transitions starting from state s1 and any ending status until this age */ - for(jk=1,pos=0, pospropta=0.; jk <=nlstate ; jk++){ - pos += pp[jk]; /* pos is the total number of transitions until this age */ - posprop[jk] += prop[jk][iage]; /* prop is the number of transitions from a live state - from jk at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */ - pospropta += prop[jk][iage]; /* prop is the number of transitions from a live state - from jk at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */ + for(s1=1,pos=0, pospropta=0.; s1 <=nlstate ; s1++){ + pos += pp[s1]; /* pos is the total number of transitions until this age */ + posprop[s1] += prop[s1][iage]; /* prop is the number of transitions from a live state + from s1 at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */ + pospropta += prop[s1][iage]; /* prop is the number of transitions from a live state + from s1 at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */ } - for(jk=1; jk <=nlstate ; jk++){ + + /* Writing ficresp */ + if(cptcoveff==0 && nj==1){ /* no covariate and first pass */ + if( iage <= iagemax){ + fprintf(ficresp," %d",iage); + } + }else if( nj==2){ + if( iage <= iagemax){ + fprintf(ficresp," %d",iage); + for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresp, " %d %d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); + } + } + for(s1=1; s1 <=nlstate ; s1++){ if(pos>=1.e-5){ if(first==1) - printf(" %d.=%.0f prev[%d]=%.1f%%",jk,pp[jk],jk,100*pp[jk]/pos); - fprintf(ficlog," %d.=%.0f prev[%d]=%.1f%%",jk,pp[jk],jk,100*pp[jk]/pos); + printf(" %d.=%.0f prev[%d]=%.1f%%",s1,pp[s1],s1,100*pp[s1]/pos); + fprintf(ficlog," %d.=%.0f prev[%d]=%.1f%%",s1,pp[s1],s1,100*pp[s1]/pos); }else{ if(first==1) - printf(" %d.=%.0f prev[%d]=NaNQ%%",jk,pp[jk],jk); - fprintf(ficlog," %d.=%.0f prev[%d]=NaNQ%%",jk,pp[jk],jk); + printf(" %d.=%.0f prev[%d]=NaNQ%%",s1,pp[s1],s1); + fprintf(ficlog," %d.=%.0f prev[%d]=NaNQ%%",s1,pp[s1],s1); } if( iage <= iagemax){ if(pos>=1.e-5){ - fprintf(ficresp," %d %.5f %.0f %.0f",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta); - /* fprintf(ficresp, "%d %d %d %.5f %.0f %.0f",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)],iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta); */ - - fprintf(ficresphtm,"",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta); - /*probs[iage][jk][j1]= pp[jk]/pos;*/ - /*printf("\niage=%d jk=%d j1=%d %.5f %.0f %.0f %f",iage,jk,j1,pp[jk]/pos, pp[jk],pos,probs[iage][jk][j1]);*/ - } - else{ - fprintf(ficresp," %d NaNq %.0f %.0f",iage,prop[jk][iage],pospropta); - fprintf(ficresphtm,"",iage, prop[jk][iage],pospropta); + if(cptcoveff==0 && nj==1){ /* no covariate and first pass */ + fprintf(ficresp," %.5f %.0f %.0f",prop[s1][iage]/pospropta, prop[s1][iage],pospropta); + }else if( nj==2){ + fprintf(ficresp," %.5f %.0f %.0f",prop[s1][iage]/pospropta, prop[s1][iage],pospropta); + } + fprintf(ficresphtm,"",iage,prop[s1][iage]/pospropta, prop[s1][iage],pospropta); + /*probs[iage][s1][j1]= pp[s1]/pos;*/ + /*printf("\niage=%d s1=%d j1=%d %.5f %.0f %.0f %f",iage,s1,j1,pp[s1]/pos, pp[s1],pos,probs[iage][s1][j1]);*/ + } else{ + if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp," NaNq %.0f %.0f",prop[s1][iage],pospropta); + fprintf(ficresphtm,"",iage, prop[s1][iage],pospropta); } } - pospropt[jk] +=posprop[jk]; - } /* end loop jk */ + pospropt[s1] +=posprop[s1]; + } /* end loop s1 */ /* pospropt=0.; */ - for(jk=-1; jk <=nlstate+ndeath; jk++){ + for(s1=-1; s1 <=nlstate+ndeath; s1++){ for(m=-1; m <=nlstate+ndeath; m++){ - if(freq[jk][m][iage] !=0 ) { /* minimizing output */ + if(freq[s1][m][iage] !=0 ) { /* minimizing output */ if(first==1){ - printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]); + printf(" %d%d=%.0f",s1,m,freq[s1][m][iage]); } - /* printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]); */ - fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iage]); + /* printf(" %d%d=%.0f",s1,m,freq[s1][m][iage]); */ + fprintf(ficlog," %d%d=%.0f",s1,m,freq[s1][m][iage]); } - if(jk!=0 && m!=0) - fprintf(ficresphtmfr," ",freq[jk][m][iage]); + if(s1!=0 && m!=0) + fprintf(ficresphtmfr," ",freq[s1][m][iage]); } - } /* end loop jk */ + } /* end loop s1 */ posproptt=0.; - for(jk=1; jk <=nlstate; jk++){ - posproptt += pospropt[jk]; + for(s1=1; s1 <=nlstate; s1++){ + posproptt += pospropt[s1]; } fprintf(ficresphtmfr,"\n "); - if(iage <= iagemax){ - fprintf(ficresp,"\n"); - fprintf(ficresphtm,"\n"); + fprintf(ficresphtm,"\n"); + if((cptcoveff==0 && nj==1)|| nj==2 ) { + if(iage <= iagemax) + fprintf(ficresp,"\n"); } if(first==1) printf("Others in log...\n"); fprintf(ficlog,"\n"); } /* end loop age iage */ + fprintf(ficresphtm,""); - for(jk=1; jk <=nlstate ; jk++){ + for(s1=1; s1 <=nlstate ; s1++){ if(posproptt < 1.e-5){ - fprintf(ficresphtm,"",pospropt[jk],posproptt); + fprintf(ficresphtm,"",pospropt[s1],posproptt); }else{ - fprintf(ficresphtm,"",pospropt[jk]/posproptt,pospropt[jk],posproptt); + fprintf(ficresphtm,"",pospropt[s1]/posproptt,pospropt[s1],posproptt); } } fprintf(ficresphtm,"\n"); @@ -4687,66 +4735,66 @@ Title=%s
Datafile=%s Firstpass=%d La fprintf(ficlog,"\n"); if(j!=0){ printf("#Freqsummary: Starting values for combination j1=%d:\n", j1); - for(i=1,jk=1; i <=nlstate; i++){ + for(i=1,s1=1; i <=nlstate; i++){ for(k=1; k <=(nlstate+ndeath); k++){ if (k != i) { - for(jj=1; jj <=ncovmodel; jj++){ /* For counting jk */ + for(jj=1; jj <=ncovmodel; jj++){ /* For counting s1 */ if(jj==1){ /* Constant case (in fact cste + age) */ if(j1==1){ /* All dummy covariates to zero */ freq[i][k][iagemax+4]=freq[i][k][iagemax+3]; /* Stores case 0 0 0 */ freq[i][i][iagemax+4]=freq[i][i][iagemax+3]; /* Stores case 0 0 0 */ printf("%d%d ",i,k); fprintf(ficlog,"%d%d ",i,k); - printf("%12.7f ln(%.0f/%.0f)= %f, OR=%f sd=%f \n",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]),freq[i][k][iagemax+3]/freq[i][i][iagemax+3], sqrt(1/freq[i][k][iagemax+3]+1/freq[i][i][iagemax+3])); - fprintf(ficlog,"%12.7f ln(%.0f/%.0f)= %12.7f \n",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])); - pstart[jk]= log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]); + printf("%12.7f ln(%.0f/%.0f)= %f, OR=%f sd=%f \n",p[s1],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]),freq[i][k][iagemax+3]/freq[i][i][iagemax+3], sqrt(1/freq[i][k][iagemax+3]+1/freq[i][i][iagemax+3])); + fprintf(ficlog,"%12.7f ln(%.0f/%.0f)= %12.7f \n",p[s1],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])); + pstart[s1]= log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]); } }else if((j1==1) && (jj==2 || nagesqr==1)){ /* age or age*age parameter without covariate V4*age (to be done later) */ for(iage=iagemin; iage <= iagemax+3; iage++){ x[iage]= (double)iage; y[iage]= log(freq[i][k][iage]/freq[i][i][iage]); - /* printf("i=%d, k=%d, jk=%d, j1=%d, jj=%d, y[%d]=%f\n",i,k,jk,j1,jj, iage, y[iage]); */ + /* printf("i=%d, k=%d, s1=%d, j1=%d, jj=%d, y[%d]=%f\n",i,k,s1,j1,jj, iage, y[iage]); */ } linreg(iagemin,iagemax,&no,x,y,&a,&b,&r, &sa, &sb ); /* y= a+b*x with standard errors */ - pstart[jk]=b; - pstart[jk-1]=a; + pstart[s1]=b; + pstart[s1-1]=a; }else if( j1!=1 && (j1==2 || (log(j1-1.)/log(2.)-(int)(log(j1-1.)/log(2.))) <0.010) && ( TvarsDind[(int)(log(j1-1.)/log(2.))+1]+2+nagesqr == jj) && Dummy[jj-2-nagesqr]==0){ /* We want only if the position, jj, in model corresponds to unique covariate equal to 1 in j1 combination */ printf("j1=%d, jj=%d, (int)(log(j1-1.)/log(2.))+1=%d, TvarsDind[(int)(log(j1-1.)/log(2.))+1]=%d\n",j1, jj,(int)(log(j1-1.)/log(2.))+1,TvarsDind[(int)(log(j1-1.)/log(2.))+1]); printf("j1=%d, jj=%d, (log(j1-1.)/log(2.))+1=%f, TvarsDind[(int)(log(j1-1.)/log(2.))+1]=%d\n",j1, jj,(log(j1-1.)/log(2.))+1,TvarsDind[(int)(log(j1-1.)/log(2.))+1]); - pstart[jk]= log((freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4])); + pstart[s1]= log((freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4])); printf("%d%d ",i,k); fprintf(ficlog,"%d%d ",i,k); - printf("jk=%d,i=%d,k=%d,p[%d]=%12.7f ln((%.0f/%.0f)/(%.0f/%.0f))= %f, OR=%f sd=%f \n",jk,i,k,jk,p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3],freq[i][k][iagemax+4],freq[i][i][iagemax+4], log((freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4])),(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4]), sqrt(1/freq[i][k][iagemax+3]+1/freq[i][i][iagemax+3]+1/freq[i][k][iagemax+4]+1/freq[i][i][iagemax+4])); + printf("s1=%d,i=%d,k=%d,p[%d]=%12.7f ln((%.0f/%.0f)/(%.0f/%.0f))= %f, OR=%f sd=%f \n",s1,i,k,s1,p[s1],freq[i][k][iagemax+3],freq[i][i][iagemax+3],freq[i][k][iagemax+4],freq[i][i][iagemax+4], log((freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4])),(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4]), sqrt(1/freq[i][k][iagemax+3]+1/freq[i][i][iagemax+3]+1/freq[i][k][iagemax+4]+1/freq[i][i][iagemax+4])); }else{ /* Other cases, like quantitative fixed or varying covariates */ ; } /* printf("%12.7f )", param[i][jj][k]); */ /* fprintf(ficlog,"%12.7f )", param[i][jj][k]); */ - jk++; + s1++; } /* end jj */ } /* end k!= i */ } /* end k */ - } /* end i, jk */ + } /* end i, s1 */ } /* end j !=0 */ } /* end selected combination of covariate j1 */ if(j==0){ /* We can estimate starting values from the occurences in each case */ printf("#Freqsummary: Starting values for the constants:\n"); fprintf(ficlog,"\n"); - for(i=1,jk=1; i <=nlstate; i++){ + for(i=1,s1=1; i <=nlstate; i++){ for(k=1; k <=(nlstate+ndeath); k++){ if (k != i) { printf("%d%d ",i,k); fprintf(ficlog,"%d%d ",i,k); for(jj=1; jj <=ncovmodel; jj++){ - pstart[jk]=p[jk]; /* Setting pstart to p values by default */ + pstart[s1]=p[s1]; /* Setting pstart to p values by default */ if(jj==1){ /* Age has to be done */ - pstart[jk]= log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]); - printf("%12.7f ln(%.0f/%.0f)= %12.7f ",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])); - fprintf(ficlog,"%12.7f ln(%.0f/%.0f)= %12.7f ",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])); + pstart[s1]= log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]); + printf("%12.7f ln(%.0f/%.0f)= %12.7f ",p[s1],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])); + fprintf(ficlog,"%12.7f ln(%.0f/%.0f)= %12.7f ",p[s1],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])); } /* printf("%12.7f )", param[i][jj][k]); */ /* fprintf(ficlog,"%12.7f )", param[i][jj][k]); */ - jk++; + s1++; } printf("\n"); fprintf(ficlog,"\n"); @@ -4755,17 +4803,17 @@ Title=%s
Datafile=%s Firstpass=%d La } printf("#Freqsummary\n"); fprintf(ficlog,"\n"); - for(jk=-1; jk <=nlstate+ndeath; jk++){ - for(m=-1; m <=nlstate+ndeath; m++){ - /* param[i]|j][k]= freq[jk][m][iagemax+3] */ - printf(" %d%d=%.0f",jk,m,freq[jk][m][iagemax+3]); - fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iagemax+3]); - /* if(freq[jk][m][iage] !=0 ) { /\* minimizing output *\/ */ - /* printf(" %d%d=%.0f",jk,m,freq[jk][m][iagemax+3]); */ - /* fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iagemax+3]); */ + for(s1=-1; s1 <=nlstate+ndeath; s1++){ + for(s2=-1; s2 <=nlstate+ndeath; s2++){ + /* param[i]|j][k]= freq[s1][s2][iagemax+3] */ + printf(" %d%d=%.0f",s1,s2,freq[s1][s2][iagemax+3]); + fprintf(ficlog," %d%d=%.0f",s1,s2,freq[s1][s2][iagemax+3]); + /* if(freq[s1][s2][iage] !=0 ) { /\* minimizing output *\/ */ + /* printf(" %d%d=%.0f",s1,s2,freq[s1][s2][iagemax+3]); */ + /* fprintf(ficlog," %d%d=%.0f",s1,s2,freq[s1][s2][iagemax+3]); */ /* } */ } - } /* end loop jk */ + } /* end loop s1 */ printf("\n"); fprintf(ficlog,"\n"); @@ -6773,7 +6821,34 @@ void printinggnuplot(char fileresu[], ch if (i==cpt) fprintf(ficgp," %%lf (%%lf)"); else fprintf(ficgp," %%*lf (%%*lf)"); } - fprintf(ficgp,"\" t\"\" w l lt 1,\"%s\" every :::%d::%d u 1:($%d) t\"Observed prevalence\" w l lt 2",subdirf2(fileresu,"P_"),k1-1,k1-1,2+4*(cpt-1)); + /* fprintf(ficgp,"\" t\"\" w l lt 1,\"%s\" every :::%d::%d u 1:($%d) t\"Observed prevalence\" w l lt 2",subdirf2(fileresu,"P_"),k1-1,k1-1,2+4*(cpt-1)); */ + + fprintf(ficgp,"\" t\"\" w l lt 1,\"%s\" u 1:((",subdirf2(fileresu,"P_")); + if(cptcoveff ==0){ + fprintf(ficgp,"$%d)) t 'Observed prevalence in state %d' with line lt 3", 2+(cpt-1), cpt ); + }else{ + kl=0; + for (k=1; k<=cptcoveff; k++){ /* For each combination of covariate */ + lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate value corresponding to k1 combination and kth covariate */ + /* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */ + /* decodtabm(1,2,4) = 1 because h=1 k= 1 (1) 1 1 */ + /* decodtabm(13,3,4)= 2 because h=13 k= 1 1 (2) 2 */ + vlv= nbcode[Tvaraff[k]][lv]; + kl++; + /* kl=6+(cpt-1)*(nlstate+1)+1+(i-1); /\* 6+(1-1)*(2+1)+1+(1-1)=7, 6+(2-1)(2+1)+1+(1-1)=10 *\/ */ + /*6+(cpt-1)*(nlstate+1)+1+(i-1)+(nlstate+1)*nlstate; 6+(1-1)*(2+1)+1+(1-1) +(2+1)*2=13 */ + /*6+1+(i-1)+(nlstate+1)*nlstate; 6+1+(1-1) +(2+1)*2=13 */ + /* '' u 6:(($1==1 && $2==0 && $3==2 && $4==0)? $9/(1.-$15) : 1/0):($5==2000? 3:2) t 'p.1' with line lc variable*/ + if(k==cptcoveff){ + fprintf(ficgp,"$%d==%d && $%d==%d)? $%d : 1/0) t 'Observed prevalence in state %d' w l lt 2",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv], \ + 2+cptcoveff*2+3*(cpt-1), cpt ); /* 4 or 6 ?*/ + }else{ + fprintf(ficgp,"$%d==%d && $%d==%d && ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv]); + kl++; + } + } /* end covariate */ + } /* end if no covariate */ + if(backcast==1){ /* We need to get the corresponding values of the covariates involved in this combination k1 */ /* fprintf(ficgp,",\"%s\" every :::%d::%d u 1:($%d) t\"Backward stable prevalence\" w l lt 3",subdirf2(fileresu,"PLB_"),k1-1,k1-1,1+cpt); */ fprintf(ficgp,",\"%s\" u 1:((",subdirf2(fileresu,"PLB_")); /* Age is in 1, nres in 2 to be fixed */
Age%d%d%d%d
%d%d%.5f%.0f%.0f%dNaNq%.0f%.0f%d%.5f%.0f%.0f%dNaNq%.0f%.0f%.0f%.0f
TotNanq%.0f%.0fNanq%.0f%.0f%.5f%.0f%.0f%.5f%.0f%.0f