--- imach/src/imach.c 2022/09/11 07:58:42 1.341 +++ imach/src/imach.c 2022/09/20 00:01:38 1.348 @@ -1,6 +1,45 @@ -/* $Id: imach.c,v 1.341 2022/09/11 07:58:42 brouard Exp $ +/* $Id: imach.c,v 1.348 2022/09/20 00:01:38 brouard Exp $ $State: Exp $ $Log: imach.c,v $ + Revision 1.348 2022/09/20 00:01:38 brouard + Summary: version 0.99r43 + + * imach.c (Module): Version 0.99r42 needed a newer version of + Gnuplot. But newer version 0.99r43 should run with the Gnuplot + version 5.0 or 5.1 distributed with IMaCh. + + Revision 1.347 2022/09/18 14:36:44 brouard + Summary: version 0.99r42 + + Revision 1.346 2022/09/16 13:52:36 brouard + * src/imach.c (Module): 0.99r41 Was an error when product of timevarying and fixed. Using FixedV[of name] now. Thank you Feinuo + + Revision 1.345 2022/09/16 13:40:11 brouard + Summary: Version 0.99r41 + + * imach.c (Module): 0.99r41 Was an error when product of timevarying and fixed. Using FixedV[of name] now. Thank you Feinuo + + Revision 1.344 2022/09/14 19:33:30 brouard + Summary: version 0.99r40 + + * imach.c (Module): Fixing names of variables in T_ (thanks to Feinuo) + + Revision 1.343 2022/09/14 14:22:16 brouard + Summary: version 0.99r39 + + * imach.c (Module): Version 0.99r39 with colored dummy covariates + (fixed or time varying), using new last columns of + ILK_parameter.txt file. + + Revision 1.342 2022/09/11 19:54:09 brouard + Summary: 0.99r38 + + * imach.c (Module): Adding timevarying products of any kinds, + should work before shifting cotvar from ncovcol+nqv columns in + order to have a correspondance between the column of cotvar and + the id of column. + (Module): Some cleaning and adding covariates in ILK.txt + Revision 1.341 2022/09/11 07:58:42 brouard Summary: Version 0.99r38 @@ -1285,6 +1324,8 @@ typedef struct { #define MAXLINE 2048 /* Was 256 and 1024. Overflow with 312 with 2 states and 4 covariates. Should be ok */ #define GNUPLOTPROGRAM "gnuplot" +#define GNUPLOTVERSION 5.1 +double gnuplotversion=GNUPLOTVERSION; /*#define GNUPLOTPROGRAM "..\\gp37mgw\\wgnuplot"*/ #define FILENAMELENGTH 256 @@ -1321,15 +1362,16 @@ typedef struct { #define ODIRSEPARATOR '\\' #endif -/* $Id: imach.c,v 1.341 2022/09/11 07:58:42 brouard Exp $ */ +/* $Id: imach.c,v 1.348 2022/09/20 00:01:38 brouard Exp $ */ /* $State: Exp $ */ #include "version.h" char version[]=__IMACH_VERSION__; char copyright[]="September 2022,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2020, Nihon University 2021-202, INED 2000-2022"; -char fullversion[]="$Revision: 1.341 $ $Date: 2022/09/11 07:58:42 $"; +char fullversion[]="$Revision: 1.348 $ $Date: 2022/09/20 00:01:38 $"; char strstart[80]; char optionfilext[10], optionfilefiname[FILENAMELENGTH]; int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */ +int debugILK=0; /* debugILK is set by a #d in a comment line */ int nagesqr=0, nforce=0; /* nagesqr=1 if model is including age*age, number of forces */ /* Number of covariates model (1)=V2+V1+ V3*age+V2*V4 */ /* Model(2) V1 + V2 + V3 + V8 + V7*V8 + V5*V6 + V8*age + V3*age + age*age */ @@ -1533,8 +1575,8 @@ int **nbcode, *Tvar; /**< model=V2 => Tv # States 1=Coresidence, 2 Living alone, 3 Institution # V1=sex, V2=raedyrs Quant Fixed, State=livarnb4..livarnb11, V3=iadl4..iald11, V4=adlw4..adlw11, V5=r4bmi..r11bmi */ -/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ -/* k 1 2 3 4 5 6 7 8 9 */ +/* V5+V4+ V3+V4*V3 +V5*age+V2 +V1*V2+V1*age+V1 */ +/* kmodel 1 2 3 4 5 6 7 8 9 */ /*Typevar[k]= 0 0 0 2 1 0 2 1 0 *//*0 for simple covariate (dummy, quantitative,*/ /* fixed or varying), 1 for age product, 2 for*/ /* product */ @@ -3925,7 +3967,7 @@ double func( double *x) */ ioffset=2+nagesqr ; /* Fixed */ - for (kf=1; kf<=ncovf;kf++){ /* For each fixed covariate dummu or quant or prod */ + for (kf=1; kf<=ncovf;kf++){ /* For each fixed covariate dummy or quant or prod */ /* # V1=sex, V2=raedyrs Quant Fixed, State=livarnb4..livarnb11, V3=iadl4..iald11, V4=adlw4..adlw11, V5=r4bmi..r11bmi */ /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* TvarF[1]=Tvar[6]=2, TvarF[2]=Tvar[7]=7, TvarF[3]=Tvar[9]=1 ID of fixed covariates or product V2, V1*V2, V1 */ @@ -3954,10 +3996,10 @@ double func( double *x) for(ncovv=1, iposold=0; ncovv <= ncovvt ; ncovv++){ /* Varying covariates (single and product but no age )*/ itv=TvarVV[ncovv]; /* TvarVV={3, 1, 3} gives the name of each varying covariate */ ipos=TvarVVind[ncovv]; /* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate*/ - if(TvarFind[itv]==0){ /* Not a fixed covariate */ + if(FixedV[itv]!=0){ /* Not a fixed covariate */ cotvarv=cotvar[mw[mi][i]][TvarVV[ncovv]][i]; /* cotvar[wav][ncovcol+nqv+iv][i] */ }else{ /* fixed covariate */ - cotvarv=covar[Tvar[TvarFind[itv]]][i]; + cotvarv=covar[itv][i]; /* Error: TvarFind gives the name, that is the true column of fixed covariates, but Tvar of the model */ } if(ipos!=iposold){ /* Not a product or first of a product */ cotvarvold=cotvarv; @@ -4217,7 +4259,7 @@ double func( double *x) ipmx +=1; sw += weight[i]; 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]); */ + /* printf("num[i]=%09ld, 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",num[i],i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); */ } /* end of wave */ } /* end of individual */ }else{ /* ml=5 no inter-extrapolation no jackson =0.8a */ @@ -4356,18 +4398,40 @@ double funcone( double *x) * k= 1 2 3 4 5 6 7 8 9 * varying 1 2 3 4 5 * ncovv 1 2 3 4 5 6 7 8 - * TvarVV V3 5 1 3 3 5 1 5 + * TvarVV[ncovv] V3 5 1 3 3 5 1 5 * TvarVVind 2 3 7 7 8 8 9 9 * TvarFind[k] 1 0 0 0 0 0 0 0 0 - * cotvar starts at ntv=2 (because of V3 V4) */ + /* Other model ncovcol=5 nqv=0 ntv=3 nqtv=0 nlstate=3 + * V2 V3 V4 are fixed V6 V7 are timevarying so V8 and V5 are not in the model and product column will start at 9 Tvar[4]=6 + * FixedV[ncovcol+qv+ntv+nqtv] V5 + * V1 V2 V3 V4 V5 V6 V7 V8 + * 0 0 0 0 0 1 1 1 + * model= V2 + V3 + V4 + V6 + V7 + V6*V2 + V7*V2 + V6*V3 + V7*V3 + V6*V4 + V7*V4 + * kmodel 1 2 3 4 5 6 7 8 9 10 11 + * ncovf 1 2 3 + * ncovvt=14 1 2 3 4 5 6 7 8 9 10 11 12 13 14 + * TvarVV[1]@14 = itv {6, 7, 6, 2, 7, 2, 6, 3, 7, 3, 6, 4, 7, 4} + * TvarVVind[1]@14= {4, 5, 6, 6, 7, 7, 8, 8, 9, 9, 10, 10, 11, 11} + * TvarFind[1]@14= {1, 2, 3, 0 } + * Tvar[1]@20= {2, 3, 4, 6, 7, 9, 10, 11, 12, 13, 14} + * TvarFind[itv] 0 0 0 + * FixedV[itv] 1 1 1 0 1 0 1 0 1 0 0 + * Tvar[TvarFind[ncovf]]=[1]=2 [2]=3 [4]=4 + * Tvar[TvarFind[itv]] [0]=? ?ncovv 1 à ncovvt] + * Not a fixed cotvar[mw][itv][i] 6 7 6 2 7, 2, 6, 3, 7, 3, 6, 4, 7, 4} + * fixed covar[itv] [6] [7] [6][2] + */ + for(ncovv=1, iposold=0; ncovv <= ncovvt ; ncovv++){ /* Varying covariates (single and product but no age) including individual from products */ - itv=TvarVV[ncovv]; /* TvarVV={3, 1, 3} gives the name of each varying covariate */ + itv=TvarVV[ncovv]; /* TvarVV={3, 1, 3} gives the name of each varying covariate, exploding product */ ipos=TvarVVind[ncovv]; /* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate*/ - if(TvarFind[itv]==0){ /* Not a fixed covariate */ - cotvarv=cotvar[mw[mi][i]][TvarVV[ncovv]][i]; /* because cotvar starts now at first ncovcol+nqv+ (1 to nqtv) */ + /* if(TvarFind[itv]==0){ /\* Not a fixed covariate? Could be a fixed covariate of a product with a higher than ncovcol+nqv, itv *\/ */ + if(FixedV[itv]!=0){ /* Not a fixed covariate? Could be a fixed covariate of a product with a higher than ncovcol+nqv, itv */ + cotvarv=cotvar[mw[mi][i]][TvarVV[ncovv]][i]; /* because cotvar starts now at first ncovcol+nqv+ntv+nqtv (1 to nqtv) */ }else{ /* fixed covariate */ - cotvarv=covar[Tvar[TvarFind[itv]]][i]; + /* cotvarv=covar[Tvar[TvarFind[itv]]][i]; /\* Error: TvarFind gives the name, that is the true column of fixed covariates, but Tvar of the model *\/ */ + cotvarv=covar[itv][i]; /* Error: TvarFind gives the name, that is the true column of fixed covariates, but Tvar of the model */ } if(ipos!=iposold){ /* Not a product or first of a product */ cotvarvold=cotvarv; @@ -4473,43 +4537,53 @@ double funcone( double *x) ipmx +=1; sw += weight[i]; ll[s[mw[mi][i]][i]] += 2*weight[i]*lli; - printf("Funcone num[i]=%ld i=%6d V= ", num[i], i); - for (kf=1; kf<=ncovf;kf++){ /* Simple and product fixed covariates without age* products *//* Missing values are set to -1 but should be dropped */ - printf("%g",covar[Tvar[TvarFind[kf]]][i]); - } - for(ncovv=1, iposold=0; ncovv <= ncovvt ; ncovv++){ /* Varying covariates (single and product but no age) including individual from products */ - ipos=TvarVVind[ncovv]; /* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate*/ - if(ipos!=iposold){ /* Not a product or first of a product */ - printf(" %g",cov[ioffset+ipos]); - }else{ - printf("*"); - } - iposold=ipos; - } - for (kk=1; kk<=cptcovage;kk++) { - if(!FixedV[Tvar[Tage[kk]]]) - printf(" %g*age",covar[Tvar[Tage[kk]]][i]); - else - printf(" %g*age",cotvar[mw[mi][i]][Tvar[Tage[kk]]][i]);/* because cotvar starts now at first ncovcol+nqv+ (1 to nqtv) */ - } - printf(" s1=%1d s2=%1d mi=%1d mw=%1d dh=%3d prob=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); + /* Printing covariates values for each contribution for checking */ + /* printf("num[i]=%09ld, 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",num[i],i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); */ if(globpr){ fprintf(ficresilk,"%09ld %6.1f %6.1f %6d %2d %2d %2d %2d %3d %15.6f %8.4f %8.3f\ %11.6f %11.6f %11.6f ", \ num[i], agebegin, ageend, i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],weight[i]*gipmx/gsw, 2*weight[i]*lli,(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); - /* printf("%09ld %6.1f %6.1f %6d %2d %2d %2d %2d %3d %15.6f %8.4f %8.3f\ */ - /* %11.6f %11.6f %11.6f ", \ */ - /* num[i], agebegin, ageend, i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],weight[i]*gipmx/gsw, */ - /* 2*weight[i]*lli,(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); */ + /* printf("%09ld %6.1f %6.1f %6d %2d %2d %2d %2d %3d %15.6f %8.4f %8.3f\ */ + /* %11.6f %11.6f %11.6f ", \ */ + /* num[i], agebegin, ageend, i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],weight[i]*gipmx/gsw, */ + /* 2*weight[i]*lli,(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); */ for(k=1,llt=0.,l=0.; k<=nlstate; k++){ llt +=ll[k]*gipmx/gsw; fprintf(ficresilk," %10.6f",-ll[k]*gipmx/gsw); /* printf(" %10.6f",-ll[k]*gipmx/gsw); */ } - fprintf(ficresilk," %10.6f\n", -llt); + fprintf(ficresilk," %10.6f ", -llt); /* printf(" %10.6f\n", -llt); */ - } + /* if(debugILK){ /\* debugILK is set by a #d in a comment line *\/ */ + /* fprintf(ficresilk,"%09ld ", num[i]); */ /* not necessary */ + for (kf=1; kf<=ncovf;kf++){ /* Simple and product fixed covariates without age* products *//* Missing values are set to -1 but should be dropped */ + fprintf(ficresilk," %g",covar[Tvar[TvarFind[kf]]][i]); + } + for(ncovv=1, iposold=0; ncovv <= ncovvt ; ncovv++){ /* Varying covariates (single and product but no age) including individual from products */ + ipos=TvarVVind[ncovv]; /* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate*/ + if(ipos!=iposold){ /* Not a product or first of a product */ + fprintf(ficresilk," %g",cov[ioffset+ipos]); + /* printf(" %g",cov[ioffset+ipos]); */ + }else{ + fprintf(ficresilk,"*"); + /* printf("*"); */ + } + iposold=ipos; + } + for (kk=1; kk<=cptcovage;kk++) { + if(!FixedV[Tvar[Tage[kk]]]){ + fprintf(ficresilk," %g*age",covar[Tvar[Tage[kk]]][i]); + /* printf(" %g*age",covar[Tvar[Tage[kk]]][i]); */ + }else{ + fprintf(ficresilk," %g*age",cotvar[mw[mi][i]][Tvar[Tage[kk]]][i]);/* because cotvar starts now at first ncovcol+nqv+ (1 to nqtv) */ + /* printf(" %g*age",cotvar[mw[mi][i]][Tvar[Tage[kk]]][i]);/\* because cotvar starts now at first ncovcol+nqv+ (1 to nqtv) *\/ */ + } + } + /* printf("\n"); */ + /* } /\* End debugILK *\/ */ + fprintf(ficresilk,"\n"); + } /* End if globpr */ } /* end of wave */ } /* end of individual */ for(k=1,l=0.; k<=nlstate; k++) l += ll[k]; @@ -4519,7 +4593,7 @@ double funcone( double *x) gipmx=ipmx; gsw=sw; } -return -l; + return -l; } @@ -4530,8 +4604,9 @@ void likelione(FILE *ficres,double p[], the selection of individuals/waves and to check the exact contribution to the likelihood. Plotting could be done. - */ - int k; + */ + void pstamp(FILE *ficres); + int k, kf, kk, kvar, ncovv, iposold, ipos; if(*globpri !=0){ /* Just counts and sums, no printings */ strcpy(fileresilk,"ILK_"); @@ -4540,13 +4615,43 @@ void likelione(FILE *ficres,double p[], printf("Problem with resultfile: %s\n", fileresilk); fprintf(ficlog,"Problem with resultfile: %s\n", fileresilk); } + pstamp(ficresilk);fprintf(ficresilk,"# model=1+age+%s\n",model); fprintf(ficresilk, "#individual(line's_record) count ageb ageend 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 ageb agend 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); - fprintf(ficresilk," -2*gipw/gsw*weight*ll(total)\n"); - } + fprintf(ficresilk," -2*gipw/gsw*weight*ll(total) "); + + /* if(debugILK){ /\* debugILK is set by a #d in a comment line *\/ */ + for(kf=1;kf <= ncovf; kf++){ + fprintf(ficresilk,"V%d",Tvar[TvarFind[kf]]); + /* printf("V%d",Tvar[TvarFind[kf]]); */ + } + for(ncovv=1, iposold=0; ncovv <= ncovvt ; ncovv++){ + ipos=TvarVVind[ncovv]; /* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate */ + if(ipos!=iposold){ /* Not a product or first of a product */ + /* printf(" %d",ipos); */ + fprintf(ficresilk," V%d",TvarVV[ncovv]); + }else{ + /* printf("*"); */ + fprintf(ficresilk,"*"); + } + iposold=ipos; + } + for (kk=1; kk<=cptcovage;kk++) { + if(!FixedV[Tvar[Tage[kk]]]){ + /* printf(" %d*age(Fixed)",Tvar[Tage[kk]]); */ + fprintf(ficresilk," %d*age(Fixed)",Tvar[Tage[kk]]); + }else{ + fprintf(ficresilk," %d*age(Varying)",Tvar[Tage[kk]]);/* because cotvar starts now at first ncovcol+nqv+ (1 to nqtv) */ + /* printf(" %d*age(Varying)",Tvar[Tage[kk]]);/\* because cotvar starts now at first ncovcol+nqv+ (1 to nqtv) *\/ */ + } + } + /* } /\* End if debugILK *\/ */ + /* printf("\n"); */ + fprintf(ficresilk,"\n"); + } /* End glogpri */ *fretone=(*func)(p); if(*globpri !=0){ @@ -4558,16 +4663,68 @@ void likelione(FILE *ficres,double p[], 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,"\n
Equation of the model: model=1+age+%s
\n",model); - for (k=1; k<= nlstate ; k++) { - fprintf(fichtm,"
- Probability p%dj by origin %d and destination j. Dot's sizes are related to corresponding weight: %s-p%dj.png
\ -",k,k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k); - } 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_")); +\n",subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_")); fprintf(fichtm,"
- and by state of destination %s-dest.png
\ -",subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_")); +\n",subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_"),subdirf2(optionfilefiname,"ILK_")); + + for (k=1; k<= nlstate ; k++) { + fprintf(fichtm,"
- Probability p%dj by origin %d and destination j. Dot's sizes are related to corresponding weight: %s-p%dj.png
\n \ +\n",k,k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k); + for(kf=1; kf <= ncovf; kf++){ /* For each simple dummy covariate of the model */ + /* kvar=Tvar[TvarFind[kf]]; */ /* variable */ + fprintf(fichtm,"
- Probability p%dj by origin %d and destination j with colored covariate V%d. Same dot size of all points but with a different color for transitions with dummy variable V%d=1 at beginning of transition (keeping former color for V%d=0): %s-p%dj.png
\ +",k,k,Tvar[TvarFind[kf]],Tvar[TvarFind[kf]],Tvar[TvarFind[kf]],subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,Tvar[TvarFind[kf]]); + } + for(ncovv=1, iposold=0; ncovv <= ncovvt ; ncovv++){ /* Loop on the time varying extended covariates (with extension of Vn*Vm */ + ipos=TvarVVind[ncovv]; /* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate */ + kvar=TvarVV[ncovv]; /* TvarVV={3, 1, 3} gives the name of each varying covariate */ + /* printf("DebugILK fichtm ncovv=%d, kvar=TvarVV[ncovv]=V%d, ipos=TvarVVind[ncovv]=%d, Dummy[ipos]=%d, Typevar[ipos]=%d\n", ncovv,kvar,ipos,Dummy[ipos],Typevar[ipos]); */ + if(ipos!=iposold){ /* Not a product or first of a product */ + /* fprintf(ficresilk," V%d",TvarVV[ncovv]); */ + /* printf(" DebugILK fichtm ipos=%d != iposold=%d\n", ipos, iposold); */ + if(Dummy[ipos]==0 && Typevar[ipos]==0){ /* Only if dummy time varying: Dummy(0, 1=quant singor prod without age,2 dummy*age, 3quant*age) Typevar (0 single, 1=*age,2=Vn*vm) */ + fprintf(fichtm,"
- Probability p%dj by origin %d and destination j with colored time varying dummy covariate V%d. Same dot size of all points but with a different color for transitions with dummy variable V%d=1 at beginning of transition (keeping former color for V%d=0): %s-p%dj.png
\ +",k,k,kvar,kvar,kvar,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,kvar); + } /* End only for dummies time varying (single?) */ + }else{ /* Useless product */ + /* printf("*"); */ + /* fprintf(ficresilk,"*"); */ + } + iposold=ipos; + } /* For each time varying covariate */ + } /* End loop on states */ + +/* if(debugILK){ */ +/* for(kf=1; kf <= ncovf; kf++){ /\* For each simple dummy covariate of the model *\/ */ +/* /\* kvar=Tvar[TvarFind[kf]]; *\/ /\* variable *\/ */ +/* for (k=1; k<= nlstate ; k++) { */ +/* fprintf(fichtm,"
- Probability p%dj by origin %d and destination j with colored covariate V%. Same dot size of all points but with a different color for transitions with dummy variable V%d=1 at beginning of transition (keeping former color for V%d=0): %s-p%dj.png
\ */ +/* ",k,k,Tvar[TvarFind[kf]],Tvar[TvarFind[kf]],Tvar[TvarFind[kf]],subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,Tvar[TvarFind[kf]]); */ +/* } */ +/* } */ +/* for(ncovv=1, iposold=0; ncovv <= ncovvt ; ncovv++){ /\* Loop on the time varying extended covariates (with extension of Vn*Vm *\/ */ +/* ipos=TvarVVind[ncovv]; /\* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate *\/ */ +/* kvar=TvarVV[ncovv]; /\* TvarVV={3, 1, 3} gives the name of each varying covariate *\/ */ +/* /\* printf("DebugILK fichtm ncovv=%d, kvar=TvarVV[ncovv]=V%d, ipos=TvarVVind[ncovv]=%d, Dummy[ipos]=%d, Typevar[ipos]=%d\n", ncovv,kvar,ipos,Dummy[ipos],Typevar[ipos]); *\/ */ +/* if(ipos!=iposold){ /\* Not a product or first of a product *\/ */ +/* /\* fprintf(ficresilk," V%d",TvarVV[ncovv]); *\/ */ +/* /\* printf(" DebugILK fichtm ipos=%d != iposold=%d\n", ipos, iposold); *\/ */ +/* if(Dummy[ipos]==0 && Typevar[ipos]==0){ /\* Only if dummy time varying: Dummy(0, 1=quant singor prod without age,2 dummy*age, 3quant*age) Typevar (0 single, 1=*age,2=Vn*vm) *\/ */ +/* for (k=1; k<= nlstate ; k++) { */ +/* fprintf(fichtm,"
- Probability p%dj by origin %d and destination j. Same dot size of all points but with a different color for transitions with dummy variable V%d=1 at beginning of transition (keeping former color for V%d=0): %s-p%dj.png
\ */ +/* ",k,k,kvar,kvar,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,subdirf2(optionfilefiname,"ILK_"),k,kvar); */ +/* } /\* End state *\/ */ +/* } /\* End only for dummies time varying (single?) *\/ */ +/* }else{ /\* Useless product *\/ */ +/* /\* printf("*"); *\/ */ +/* /\* fprintf(ficresilk,"*"); *\/ */ +/* } */ +/* iposold=ipos; */ +/* } /\* For each time varying covariate *\/ */ +/* }/\* End debugILK *\/ */ fflush(fichtm); - } + }/* End globpri */ return; } @@ -6136,7 +6293,7 @@ void concatwav(int wav[], int **dh, int /* Loop on covariates without age and products and no quantitative variable */ for (k=1; k<=cptcovt; k++) { /* cptcovt: total number of covariates of the model (2) nbocc(+)+1 = 8 excepting constant and age and age*age */ for (j=-1; (j < maxncov); j++) Ndum[j]=0; - printf("Testing k=%d, cptcovt=%d\n",k, cptcovt); + /* printf("Testing k=%d, cptcovt=%d\n",k, cptcovt); */ if(Dummy[k]==0 && Typevar[k] !=1 && Typevar[k] != 2){ /* Dummy covariate and not age product nor fixed product */ switch(Fixed[k]) { case 0: /* Testing on fixed dummy covariate, simple or product of fixed */ @@ -6234,7 +6391,7 @@ void concatwav(int wav[], int **dh, int break; } /* end switch */ } /* end dummy test */ - if(Dummy[k]==1 && Typevar[k] !=1){ /* Quantitative covariate and not age product */ + if(Dummy[k]==1 && Typevar[k] !=1 && Fixed ==0){ /* Fixed Quantitative covariate and not age product */ for (i=1; i<=imx; i++) { /* Loop on individuals: reads the data file to get the maximum value of the modality of this covariate Vj*/ if(Tvar[k]<=0 || Tvar[k]>=NCOVMAX){ printf("Error k=%d \n",k); @@ -6700,6 +6857,7 @@ void concatwav(int wav[], int **dh, int /* fprintf(ficresprobmorprev," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ /* } */ for (j=1; j<= cptcovs; j++){ /* For each selected (single) quantitative value */ /* To be done*/ + /* fprintf(ficresprobmorprev," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); */ fprintf(ficresprobmorprev," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); } /* for(j=1;j<=cptcoveff;j++) */ @@ -7334,7 +7492,7 @@ To be simple, these graphs help to under for(nres=1;nres <=nresult; nres++){ /* For each resultline */ for(j1=1; j1<=tj;j1++){ /* For any combination of dummy covariates, fixed and varying */ - printf("Varprob TKresult[nres]=%d j1=%d, nres=%d, cptcovn=%d, cptcoveff=%d tj=%d cptcovs=%d\n", TKresult[nres], j1, nres, cptcovn, cptcoveff, tj, cptcovs); + /* printf("Varprob TKresult[nres]=%d j1=%d, nres=%d, cptcovn=%d, cptcoveff=%d tj=%d cptcovs=%d\n", TKresult[nres], j1, nres, cptcovn, cptcoveff, tj, cptcovs); */ if(tj != 1 && TKresult[nres]!= j1) continue; @@ -7350,7 +7508,7 @@ To be simple, these graphs help to under /* Including quantitative variables of the resultline to be done */ for (z1=1; z1<=cptcovs; z1++){ /* Loop on each variable of this resultline */ - printf("Varprob modelresult[%d][%d]=%d model=1+age+%s \n",nres, z1, modelresult[nres][z1], model); + /* printf("Varprob modelresult[%d][%d]=%d model=1+age+%s \n",nres, z1, modelresult[nres][z1], model); */ fprintf(ficlog,"Varprob modelresult[%d][%d]=%d model=1+age+%s \n",nres, z1, modelresult[nres][z1], model); /* fprintf(ficlog,"Varprob modelresult[%d][%d]=%d model=1+age+%s resultline[%d]=%s \n",nres, z1, modelresult[nres][z1], model, nres, resultline[nres]); */ if(Dummy[modelresult[nres][z1]]==0){/* Dummy variable of the variable in position modelresult in the model corresponding to z1 in resultline */ @@ -8009,7 +8167,7 @@ void printinggnuplot(char fileresu[], ch char dirfileres[132],optfileres[132]; char gplotcondition[132], gplotlabel[132]; - int cpt=0,k1=0,i=0,k=0,j=0,jk=0,k2=0,k3=0,k4=0,ij=0, ijp=0, l=0; + int cpt=0,k1=0,i=0,k=0,j=0,jk=0,k2=0,k3=0,k4=0,kf=0,kvar=0,kk=0,ipos=0,iposold=0,ij=0, ijp=0, l=0; int lv=0, vlv=0, kl=0; int ng=0; int vpopbased; @@ -8035,7 +8193,7 @@ void printinggnuplot(char fileresu[], ch fprintf(ficgp,"yoff=(%d > 2? 0:1);\n",nlstate); fprintf(ficgp,"\n#Peripheral arrows\nset for [i=1:%d] for [j=1:%d] arrow i*10+j from cos(pi*((1-(%d/2)*2./%d)/2+(i-1)*2./%d))-(i!=j?(i-j)/abs(i-j)*delta:0), yoff +sin(pi*((1-(%d/2)*2./%d)/2+(i-1)*2./%d)) + (i!=j?(i-j)/abs(i-j)*delta:0) rto -0.95*(cos(pi*((1-(%d/2)*2./%d)/2+(i-1)*2./%d))+(i!=j?(i-j)/abs(i-j)*delta:0) - cos(pi*((1-(%d/2)*2./%d)/2+(j-1)*2./%d)) + (i!=j?(i-j)/abs(i-j)*delta2:0)), -0.95*(sin(pi*((1-(%d/2)*2./%d)/2+(i-1)*2./%d)) + (i!=j?(i-j)/abs(i-j)*delta:0) - sin(pi*((1-(%d/2)*2./%d)/2+(j-1)*2./%d))+( i!=j?(i-j)/abs(i-j)*delta2:0)) ls (i < j? 1:2)\n",nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate); - fprintf(ficgp,"\n#Centripete arrows (turning in other direction (1-i) instead of (i-1)) \nset for [i=1:%d] arrow (%d+1)*10+i from cos(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d))-(i!=j?(i-j)/abs(i-j)*delta:0), yoff +sin(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d)) + (i!=j?(i-j)/abs(i-j)*delta:0) rto -0.80*(cos(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d))+(i!=j?(i-j)/abs(i-j)*delta:0) ), -0.80*(sin(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d)) + (i!=j?(i-j)/abs(i-j)*delta:0) + yoff ) ls 4\n",nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate); + fprintf(ficgp,"\n#Centripete arrows (turning in other direction (1-i) instead of (i-1)) \nset for [i=1:%d] for [j=1:%d] arrow (%d+1)*10+i from cos(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d))-(i!=j?(i-j)/abs(i-j)*delta:0), yoff +sin(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d)) + (i!=j?(i-j)/abs(i-j)*delta:0) rto -0.80*(cos(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d))+(i!=j?(i-j)/abs(i-j)*delta:0) ), -0.80*(sin(pi*((1-(%d/2)*2./%d)/2+(1-i)*2./%d)) + (i!=j?(i-j)/abs(i-j)*delta:0) + yoff ) ls 4\n",nlstate, nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate); fprintf(ficgp,"\n#show arrow\nunset label\n"); fprintf(ficgp,"\n#States labels, starting from 2 (2-i) instead of (1-i), was (i-1)\nset for [i=1:%d] label i sprintf(\"State %%d\",i) center at cos(pi*((1-(%d/2)*2./%d)/2+(2-i)*2./%d)), yoff+sin(pi*((1-(%d/2)*2./%d)/2+(2-i)*2./%d)) font \"helvetica, 16\" tc rgbcolor \"blue\"\n",nlstate,nlstate,nlstate,nlstate,nlstate,nlstate,nlstate); fprintf(ficgp,"\nset label %d+1 sprintf(\"State %%d\",%d+1) center at 0.,0. font \"helvetica, 16\" tc rgbcolor \"red\"\n",nlstate,nlstate); @@ -8072,6 +8230,87 @@ void printinggnuplot(char fileresu[], ch fprintf(ficgp,"\nset out;unset log\n"); /* fprintf(ficgp,"\nset out \"%s.svg\"; replot; set out; # bug gnuplot",subdirf2(optionfilefiname,"ILK_")); */ + /* Plot the probability implied in the likelihood by covariate value */ + fprintf(ficgp,"\nset ter pngcairo size 640, 480"); + /* if(debugILK==1){ */ + for(kf=1; kf <= ncovf; kf++){ /* For each simple dummy covariate of the model */ + kvar=Tvar[TvarFind[kf]]; /* variable name */ + /* k=18+Tvar[TvarFind[kf]];/\*offset because there are 18 columns in the ILK_ file but could be placed else where *\/ */ + k=18+kf;/*offset because there are 18 columns in the ILK_ file */ + for (i=1; i<= nlstate ; i ++) { + fprintf(ficgp,"\nset out \"%s-p%dj-%d.png\";set ylabel \"Probability for each individual/wave\";",subdirf2(optionfilefiname,"ILK_"),i,kvar); + fprintf(ficgp,"unset log;\n# For each simple dummy covariate of the model \n plot \"%s\"",subdirf(fileresilk)); + if(gnuplotversion >=5.2){ /* Former gnuplot versions do not have variable pointsize!! */ + fprintf(ficgp," u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? 7 : 9):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt variable ps 0.4 lc variable \\\n",i,1,k,k,i,1,kvar); + for (j=2; j<= nlstate+ndeath ; j ++) { + fprintf(ficgp,",\\\n \"\" u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? 7 : 9):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt variable ps 0.4 lc variable ",i,j,k,k,i,j,kvar); + } + }else{ + fprintf(ficgp," u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt 7 ps 0.4 lc variable \\\n",i,1,k,i,1,kvar); + for (j=2; j<= nlstate+ndeath ; j ++) { + fprintf(ficgp,",\\\n \"\" u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt 7 ps 0.4 lc variable ",i,j,k,i,j,kvar); + } + } + fprintf(ficgp,";\nset out; unset ylabel;\n"); + } + } /* End of each covariate dummy */ + for(ncovv=1, iposold=0, kk=0; ncovv <= ncovvt ; ncovv++){ + /* Other example V1 + V3 + V5 + age*V1 + age*V3 + age*V5 + V1*V3 + V3*V5 + V1*V5 + * kmodel = 1 2 3 4 5 6 7 8 9 + * varying 1 2 3 4 5 + * ncovv 1 2 3 4 5 6 7 8 + * TvarVV[ncovv] V3 5 1 3 3 5 1 5 + * TvarVVind[ncovv]=kmodel 2 3 7 7 8 8 9 9 + * TvarFind[kmodel] 1 0 0 0 0 0 0 0 0 + * kdata ncovcol=[V1 V2] nqv=0 ntv=[V3 V4] nqtv=V5 + * Dummy[kmodel] 0 0 1 2 2 3 1 1 1 + */ + ipos=TvarVVind[ncovv]; /* TvarVVind={2, 5, 5] gives the position in the model of the ncovv th varying covariate */ + kvar=TvarVV[ncovv]; /* TvarVV={3, 1, 3} gives the name of each varying covariate */ + /* printf("DebugILK ficgp ncovv=%d, kvar=TvarVV[ncovv]=%d, ipos=TvarVVind[ncovv]=%d, Dummy[ipos]=%d, Typevar[ipos]=%d\n", ncovv,kvar,ipos,Dummy[ipos],Typevar[ipos]); */ + if(ipos!=iposold){ /* Not a product or first of a product */ + /* printf(" %d",ipos); */ + /* fprintf(ficresilk," V%d",TvarVV[ncovv]); */ + /* printf(" DebugILK ficgp suite ipos=%d != iposold=%d\n", ipos, iposold); */ + kk++; /* Position of the ncovv column in ILK_ */ + k=18+ncovf+kk; /*offset because there are 18 columns in the ILK_ file plus ncovf fixed covariate */ + if(Dummy[ipos]==0 && Typevar[ipos]==0){ /* Only if dummy time varying: Dummy(0, 1=quant singor prod without age,2 dummy*age, 3quant*age) Typevar (0 single, 1=*age,2=Vn*vm) */ + for (i=1; i<= nlstate ; i ++) { + fprintf(ficgp,"\nset out \"%s-p%dj-%d.png\";set ylabel \"Probability for each individual/wave\";",subdirf2(optionfilefiname,"ILK_"),i,kvar); + fprintf(ficgp,"unset log;\n# For each simple dummy covariate of the model \n plot \"%s\"",subdirf(fileresilk)); + + /* printf("Before DebugILK gnuplotversion=%g >=5.2\n",gnuplotversion); */ + if(gnuplotversion >=5.2){ /* Former gnuplot versions do not have variable pointsize!! */ + /* printf("DebugILK gnuplotversion=%g >=5.2\n",gnuplotversion); */ + fprintf(ficgp," u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? 7 : 9):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt variable ps 0.4 lc variable \\\n",i,1,k,k,i,1,kvar); + for (j=2; j<= nlstate+ndeath ; j ++) { + fprintf(ficgp,",\\\n \"\" u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? 7 : 9):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt variable ps 0.4 lc variable ",i,j,k,k,i,j,kvar); + } + }else{ + /* printf("DebugILK gnuplotversion=%g <5.2\n",gnuplotversion); */ + fprintf(ficgp," u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt 7 ps 0.4 lc variable \\\n",i,1,k,i,1,kvar); + for (j=2; j<= nlstate+ndeath ; j ++) { + fprintf(ficgp,",\\\n \"\" u 2:($5 == %d && $6==%d ? $10 : 1/0):($%d==0 ? $6 : $6+4) t \"p%d%d V%d\" with points pt 7 ps 0.4 lc variable ",i,j,k,i,j,kvar); + } + } + fprintf(ficgp,";\nset out; unset ylabel;\n"); + } + }/* End if dummy varying */ + }else{ /*Product */ + /* printf("*"); */ + /* fprintf(ficresilk,"*"); */ + } + iposold=ipos; + } /* For each time varying covariate */ + /* } /\* debugILK==1 *\/ */ + /* unset log; plot "rrtest1_sorted_4/ILK_rrtest1_sorted_4.txt" u 2:($4 == 1 && $5==2 ? $9 : 1/0):5 t "p12" with points lc variable */ + /* fprintf(ficgp,"\nset log y;plot \"%s\" u 2:(-$11):3 t \"All sample, all transitions\" with dots lc variable",subdirf(fileresilk)); */ + /* fprintf(ficgp,"\nreplot \"%s\" u 2:($3 <= 3 ? -$11 : 1/0):3 t \"First 3 individuals\" with line lc variable", subdirf(fileresilk)); */ + fprintf(ficgp,"\nset out;unset log\n"); + /* fprintf(ficgp,"\nset out \"%s.svg\"; replot; set out; # bug gnuplot",subdirf2(optionfilefiname,"ILK_")); */ + + + strcpy(dirfileres,optionfilefiname); strcpy(optfileres,"vpl"); /* 1eme*/ @@ -10430,9 +10669,13 @@ int readdata(char datafile[], int firsto errno=0; lval=strtol(strb,&endptr,10); /* if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN))*/ - if( strb[0]=='\0' || (*endptr != '\0')){ - printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,j,maxwav); - fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,j,maxwav);fflush(ficlog); + if( strb[0]=='\0' || (*endptr != '\0' )){ + printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,j,maxwav); + fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,j,maxwav);fflush(ficlog); + return 1; + }else if( lval==0 || lval > nlstate+ndeath){ + printf("Error in data around '%s' at line number %d for individual %d, '%s'\n Should be a state at wave %d. A state should be 1 to %d and not %ld.\n Fix your data file '%s'! Exiting.\n", strb, linei,i,line,j,nlstate+ndeath, lval, datafile);fflush(stdout); + fprintf(ficlog,"Error in data around '%s' at line number %d for individual %d, '%s'\n Should be a state at wave %d. A state should be 1 to %d and not %ld.\n Fix your data file '%s'! Exiting.\n", strb, linei,i,line,j,nlstate+ndeath, lval, datafile); fflush(ficlog); return 1; } } @@ -10631,7 +10874,7 @@ int decoderesult( char resultline[], int printf("decoderesult:%s\n",resultline); strcpy(resultsav,resultline); - printf("Decoderesult resultsav=\"%s\" resultline=\"%s\"\n", resultsav, resultline); + /* printf("Decoderesult resultsav=\"%s\" resultline=\"%s\"\n", resultsav, resultline); */ if (strlen(resultsav) >1){ j=nbocc(resultsav,'='); /**< j=Number of covariate values'=' in this resultline */ } @@ -10691,7 +10934,7 @@ int decoderesult( char resultline[], int if(Tvar[k1]==Tvarsel[k2]) {/* Tvar is coming from the model, Tvarsel from the result. Tvar[1]=5 == Tvarsel[2]=5 */ modelresult[nres][k2]=k1;/* we found a Vn=1 corrresponding to Vn*age in the model modelresult[2]=1 modelresult[1]=2 modelresult[3]=3 modelresult[6]=4 modelresult[9]=5 */ resultmodel[nres][k1]=k2; /* Added here */ - printf("Decoderesult first modelresult[k2=%d]=%d (k1) V%d*AGE\n",k2,k1,Tvar[k1]); + /* printf("Decoderesult first modelresult[k2=%d]=%d (k1) V%d*AGE\n",k2,k1,Tvar[k1]); */ match=1; /* modelresult of k2 variable of resultline is identical to k1 variable of the model good */ break; } @@ -10704,11 +10947,11 @@ int decoderesult( char resultline[], int }else if(Typevar[k1]==2){ /* Product No age We want to get the position in the resultline of the product in the model line*/ /* resultmodel[nres][of such a Vn * Vm product k1] is not unique, so can't exist, we feed Tvard[k1][1] and [2] */ match=0; - printf("Decoderesult very first Product Tvardk[k1=%d][1]=%d Tvardk[k1=%d][2]=%d V%d * V%d\n",k1,Tvardk[k1][1],k1,Tvardk[k1][2],Tvardk[k1][1],Tvardk[k1][2]); + /* printf("Decoderesult very first Product Tvardk[k1=%d][1]=%d Tvardk[k1=%d][2]=%d V%d * V%d\n",k1,Tvardk[k1][1],k1,Tvardk[k1][2],Tvardk[k1][1],Tvardk[k1][2]); */ for(k2=1; k2 <=j;k2++){/* Loop on resultline. In result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */ if(Tvardk[k1][1]==Tvarsel[k2]) {/* Tvardk is coming from the model, Tvarsel from the result. Tvar[1]=5 == Tvarsel[2]=5 */ /* modelresult[k2]=k1; */ - printf("Decoderesult first Product modelresult[k2=%d]=%d (k1) V%d * \n",k2,k1,Tvarsel[k2]); + /* printf("Decoderesult first Product modelresult[k2=%d]=%d (k1) V%d * \n",k2,k1,Tvarsel[k2]); */ match=1; /* modelresult of k2 variable of resultline is identical to k1 variable of the model good */ } } @@ -10721,7 +10964,7 @@ int decoderesult( char resultline[], int for(k2=1; k2 <=j;k2++){/* Loop on resultline. In result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */ if(Tvardk[k1][2]==Tvarsel[k2]) {/* Tvardk is coming from the model, Tvarsel from the result. Tvar[1]=5 == Tvarsel[2]=5 */ /* modelresult[k2]=k1;*/ - printf("Decoderesult second Product modelresult[k2=%d]=%d (k1) * V%d \n ",k2,k1,Tvarsel[k2]); + /* printf("Decoderesult second Product modelresult[k2=%d]=%d (k1) * V%d \n ",k2,k1,Tvarsel[k2]); */ match=1; /* modelresult of k2 variable of resultline is identical to k1 variable of the model good */ break; } @@ -10781,7 +11024,7 @@ int decoderesult( char resultline[], int /* k counting number of combination of single dummies in the equation model */ /* k4 counting single dummies in the equation model */ /* k4q counting single quantitatives in the equation model */ - if( Dummy[k1]==0 && Typevar[k1]==0 ){ /* Dummy and Single, k1 is sorting according to MODEL, but k3 to resultline */ + if( Dummy[k1]==0 && Typevar[k1]==0 ){ /* Dummy and Single, fixed or timevarying, k1 is sorting according to MODEL, but k3 to resultline */ /* k4+1= (not always if quant in model) position in the resultline V(Tvarsel)=Tvalsel=Tresult[nres][pos](value); V(Tvresult[nres][pos] (variable): V(variable)=value) */ /* modelresult[k3]=k1: k3th position in the result line corresponds to the k1 position in the model line (doesn't work with products)*/ /* Value in the (current nres) resultline of the variable at the k1th position in the model equation resultmodel[nres][k1]= k3 */ @@ -10805,7 +11048,7 @@ int decoderesult( char resultline[], int Tvresult[nres][k3]=(int)Tvarsel[k3];/* Tvresult[nres][1]=4 Tvresult[nres][3]=1 */ Tinvresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* Tinvresult[nres][4]=1 */ precov[nres][k1]=Tvalsel[k3]; /* Value from resultline of the variable at the k1 position in the model */ - printf("Decoderesult Dummy k=%d, k1=%d precov[nres=%d][k1=%d]=%.f V(k2=V%d)= Tvalsel[%d]=%d, 2**(%d)\n",k, k1, nres, k1,precov[nres][k1], k2, k3, (int)Tvalsel[k3], k4); + /* printf("Decoderesult Dummy k=%d, k1=%d precov[nres=%d][k1=%d]=%.f V(k2=V%d)= Tvalsel[%d]=%d, 2**(%d)\n",k, k1, nres, k1,precov[nres][k1], k2, k3, (int)Tvalsel[k3], k4); */ k4++;; }else if( Dummy[k1]==1 && Typevar[k1]==0 ){ /* Quantitative and single */ /* Tqresult[nres][result_position]= value of the variable at the result_position in the nres resultline */ @@ -10823,7 +11066,7 @@ int decoderesult( char resultline[], int Tqinvresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* Tqinvresult[nres][5]=25.1 */ TinvDoQresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* Tqinvresult[nres][5]=25.1 */ precov[nres][k1]=Tvalsel[k3q]; - printf("Decoderesult Quantitative nres=%d,precov[nres=%d][k1=%d]=%.f V(k2q=V%d)= Tvalsel[%d]=%d, Tvarsel[%d]=%f\n",nres, nres, k1,precov[nres][k1], k2q, k3q, Tvarsel[k3q], k3q, Tvalsel[k3q]); + /* printf("Decoderesult Quantitative nres=%d,precov[nres=%d][k1=%d]=%.f V(k2q=V%d)= Tvalsel[%d]=%d, Tvarsel[%d]=%f\n",nres, nres, k1,precov[nres][k1], k2q, k3q, Tvarsel[k3q], k3q, Tvalsel[k3q]); */ k4q++;; }else if( Dummy[k1]==2 ){ /* For dummy with age product */ /* Tvar[k1]; */ /* Age variable */ @@ -10833,16 +11076,16 @@ int decoderesult( char resultline[], int k2=(int)Tvarsel[k3]; /* nres=1 k1=2=>k3=1 Tvarsel[resultmodel[2]]= Tvarsel[1] = 4=k2 (V4); k1=3=>k3=2 Tvarsel[2]=3 (V3)*/ TinvDoQresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* TinvDoQresult[nres][4]=1 */ precov[nres][k1]=Tvalsel[k3]; - printf("Decoderesult Dummy with age k=%d, k1=%d precov[nres=%d][k1=%d]=%.f Tvar[%d]=V%d k2=Tvarsel[%d]=%d Tvalsel[%d]=%d\n",k, k1, nres, k1,precov[nres][k1], k1, Tvar[k1], k3,(int)Tvarsel[k3], k3, (int)Tvalsel[k3]); + /* printf("Decoderesult Dummy with age k=%d, k1=%d precov[nres=%d][k1=%d]=%.f Tvar[%d]=V%d k2=Tvarsel[%d]=%d Tvalsel[%d]=%d\n",k, k1, nres, k1,precov[nres][k1], k1, Tvar[k1], k3,(int)Tvarsel[k3], k3, (int)Tvalsel[k3]); */ }else if( Dummy[k1]==3 ){ /* For quant with age product */ k3q= resultmodel[nres][k1]; /* resultmodel[1(V5)] = 25.1=k3q */ k2q=(int)Tvarsel[k3q]; /* Tvarsel[resultmodel[1]]= Tvarsel[1] = 4=k2 */ TinvDoQresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* TinvDoQresult[nres][5]=25.1 */ precov[nres][k1]=Tvalsel[k3q]; - printf("Decoderesult Quantitative with age nres=%d, k1=%d, precov[nres=%d][k1=%d]=%f Tvar[%d]=V%d V(k2q=%d)= Tvarsel[%d]=%d, Tvalsel[%d]=%f\n",nres, k1, nres, k1,precov[nres][k1], k1, Tvar[k1], k2q, k3q, Tvarsel[k3q], k3q, Tvalsel[k3q]); + /* printf("Decoderesult Quantitative with age nres=%d, k1=%d, precov[nres=%d][k1=%d]=%f Tvar[%d]=V%d V(k2q=%d)= Tvarsel[%d]=%d, Tvalsel[%d]=%f\n",nres, k1, nres, k1,precov[nres][k1], k1, Tvar[k1], k2q, k3q, Tvarsel[k3q], k3q, Tvalsel[k3q]); */ }else if(Typevar[k1]==2 ){ /* For product quant or dummy (not with age) */ precov[nres][k1]=TinvDoQresult[nres][Tvardk[k1][1]] * TinvDoQresult[nres][Tvardk[k1][2]]; - printf("Decoderesult Quantitative or Dummy (not with age) nres=%d k1=%d precov[nres=%d][k1=%d]=%.f V%d(=%.f) * V%d(=%.f) \n",nres, k1, nres, k1,precov[nres][k1], Tvardk[k1][1], TinvDoQresult[nres][Tvardk[k1][1]], Tvardk[k1][2], TinvDoQresult[nres][Tvardk[k1][2]]); + /* printf("Decoderesult Quantitative or Dummy (not with age) nres=%d k1=%d precov[nres=%d][k1=%d]=%.f V%d(=%.f) * V%d(=%.f) \n",nres, k1, nres, k1,precov[nres][k1], Tvardk[k1][1], TinvDoQresult[nres][Tvardk[k1][1]], Tvardk[k1][2], TinvDoQresult[nres][Tvardk[k1][2]]); */ }else{ printf("Error Decoderesult probably a product Dummy[%d]==%d && Typevar[%d]==%d\n", k1, Dummy[k1], k1, Typevar[k1]); fprintf(ficlog,"Error Decoderesult probably a product Dummy[%d]==%d && Typevar[%d]==%d\n", k1, Dummy[k1], k1, Typevar[k1]); @@ -10931,23 +11174,23 @@ int decodemodel( char model[], int lasto * Model V2 + V1 + V3*age + V3 + V5*V6 + V7*V8 + V8*age + V8 * k= 1 2 3 4 5 6 7 8 * cptcovn number of covariates (not including constant and age ) = # of + plus 1 = 7+1=8 - * covar[k,i], value of kth covariate if not including age for individual i: + * covar[k,i], are for fixed covariates, value of kth covariate if not including age for individual i: * covar[1][i]= (V1), covar[4][i]=(V4), covar[8][i]=(V8) * Tvar[k] # of the kth covariate: Tvar[1]=2 Tvar[2]=1 Tvar[4]=3 Tvar[8]=8 * if multiplied by age: V3*age Tvar[3=V3*age]=3 (V3) Tvar[7]=8 and * Tage[++cptcovage]=k - * if products, new covar are created after ncovcol with k1 + * if products, new covar are created after ncovcol + nqv (quanti fixed) with k1 * Tvar[k]=ncovcol+k1; # of the kth covariate product: Tvar[5]=ncovcol+1=10 Tvar[6]=ncovcol+1=11 * Tprod[k1]=k; Tprod[1]=5 Tprod[2]= 6; gives the position of the k1th product * Tvard[k1][1]=m Tvard[k1][2]=m; Tvard[1][1]=5 (V5) Tvard[1][2]=6 Tvard[2][1]=7 (V7) Tvard[2][2]=8 * Tvar[cptcovn+k2]=Tvard[k1][1];Tvar[cptcovn+k2+1]=Tvard[k1][2]; * Tvar[8+1]=5;Tvar[8+2]=6;Tvar[8+3]=7;Tvar[8+4]=8 inverted * V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 V11 - * < ncovcol=8 > + * < ncovcol=8 8 fixed covariate. Additional starts at 9 (V5*V6) and 10(V7*V8) > * Model V2 + V1 + V3*age + V3 + V5*V6 + V7*V8 + V8*age + V8 d1 d1 d2 d2 * k= 1 2 3 4 5 6 7 8 9 10 11 12 - * Tvar[k]= 2 1 3 3 10 11 8 8 5 6 7 8 - * p Tvar[1]@12={2, 1, 3, 3, 11, 10, 8, 8, 7, 8, 5, 6} + * Tvard[k]= 2 1 3 3 10 11 8 8 5 6 7 8 + * p Tvar[1]@12={2, 1, 3, 3, 9, 10, 8, 8} * p Tprod[1]@2={ 6, 5} *p Tvard[1][1]@4= {7, 8, 5, 6} * covar[k][i]= V2 V1 ? V3 V5*V6? V7*V8? ? V8 @@ -11091,8 +11334,9 @@ Dummy[k] 0=dummy (0 1), 1 quantitative ( Typevar: 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product \n\ Fixed[k] 0=fixed (product or simple), 1 varying, 2 fixed with age product, 3 varying with age product \n\ Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product\n",model); - for(k=-1;k<=cptcovt; k++){ Fixed[k]=0; Dummy[k]=0;} - for(k=1, ncovf=0, nsd=0, nsq=0, ncovv=0, ncova=0, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0, ncovvt=0;k<=cptcovt; k++){ /* or cptocvt */ + for(k=-1;k<=NCOVMAX; k++){ Fixed[k]=0; Dummy[k]=0;} + for(k=1;k<=NCOVMAX; k++){TvarFind[k]=0; TvarVind[k]=0;} + for(k=1, ncovf=0, nsd=0, nsq=0, ncovv=0, ncova=0, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0, ncovvt=0;k<=cptcovt; k++){ /* or cptocvt loop on k from model */ if (Tvar[k] <=ncovcol && Typevar[k]==0 ){ /* Simple fixed dummy (<=ncovcol) covariates */ Fixed[k]= 0; Dummy[k]= 0; @@ -11180,8 +11424,8 @@ Dummy[k] 0=dummy (0 1), 1 quantitative ( TvarVQind[nqtveff]=k; /* TvarVQind[1]=1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple time varying quantitative variable */ TmodelInvQind[nqtveff]=Tvar[k]- ncovcol-nqv-ntv;/* Only simple time varying quantitative variable */ /* Tmodeliqind[k]=nqtveff;/\* Only simple time varying quantitative variable *\/ */ - printf("Quasi TmodelQind[%d]=%d,Tvar[TmodelQind[%d]]=V%d, ncovcol=%d, nqv=%d, ntv=%d,Tvar[k]- ncovcol-nqv-ntv=%d\n",nqtveff,k,nqtveff,Tvar[k], ncovcol, nqv, ntv, Tvar[k]- ncovcol-nqv-ntv); - printf("Quasi TmodelInvQind[%d]=%d\n",k,Tvar[k]- ncovcol-nqv-ntv); + /* printf("Quasi TmodelQind[%d]=%d,Tvar[TmodelQind[%d]]=V%d, ncovcol=%d, nqv=%d, ntv=%d,Tvar[k]- ncovcol-nqv-ntv=%d\n",nqtveff,k,nqtveff,Tvar[k], ncovcol, nqv, ntv, Tvar[k]- ncovcol-nqv-ntv); */ + /* printf("Quasi TmodelInvQind[%d]=%d\n",k,Tvar[k]- ncovcol-nqv-ntv); */ }else if (Typevar[k] == 1) { /* product with age */ ncova++; TvarA[ncova]=Tvar[k]; @@ -11363,8 +11607,8 @@ Dummy[k] 0=dummy (0 1), 1 quantitative ( printf("Error, current version can't treat for performance reasons, Tvar[%d]=%d, Typevar[%d]=%d\n", k, Tvar[k], k, Typevar[k]); fprintf(ficlog,"Error, current version can't treat for performance reasons, Tvar[%d]=%d, Typevar[%d]=%d\n", k, Tvar[k], k, Typevar[k]); } - printf("Decodemodel, k=%d, Tvar[%d]=V%d,Typevar=%d, Fixed=%d, Dummy=%d\n",k, k,Tvar[k],Typevar[k],Fixed[k],Dummy[k]); - printf(" modell[%d].maintype=%d, modell[%d].subtype=%d\n",k,modell[k].maintype,k,modell[k].subtype); + /* printf("Decodemodel, k=%d, Tvar[%d]=V%d,Typevar=%d, Fixed=%d, Dummy=%d\n",k, k,Tvar[k],Typevar[k],Fixed[k],Dummy[k]); */ + /* printf(" modell[%d].maintype=%d, modell[%d].subtype=%d\n",k,modell[k].maintype,k,modell[k].subtype); */ fprintf(ficlog,"Decodemodel, k=%d, Tvar[%d]=V%d,Typevar=%d, Fixed=%d, Dummy=%d\n",k, k,Tvar[k],Typevar[k],Fixed[k],Dummy[k]); } /* Searching for doublons in the model */ @@ -12518,6 +12762,8 @@ int main(int argc, char *argv[]) numlinepar++; if(line[1]=='q'){ /* This #q will quit imach (the answer is q) */ z[0]=line[1]; + }else if(line[1]=='d'){ /* For debugging individual values of covariates in ficresilk */ + debugILK=1;printf("DebugILK\n"); } /* printf("****line [1] = %c \n",line[1]); */ fputs(line, stdout); @@ -13816,7 +14062,7 @@ Please run with mle=-1 to get a correct case 13: num_filled=sscanf(line,"result:%[^\n]\n",resultlineori); nresult++; /* Sum of resultlines */ - printf("Result %d: result:%s\n",nresult, resultlineori); + /* printf("Result %d: result:%s\n",nresult, resultlineori); */ /* removefirstspace(&resultlineori); */ if(strstr(resultlineori,"v") !=0){ @@ -13825,7 +14071,7 @@ Please run with mle=-1 to get a correct return 1; } trimbb(resultline, resultlineori); /* Suppressing double blank in the resultline */ - printf("Decoderesult resultline=\"%s\" resultlineori=\"%s\"\n", resultline, resultlineori); + /* printf("Decoderesult resultline=\"%s\" resultlineori=\"%s\"\n", resultline, resultlineori); */ if(nresult > MAXRESULTLINESPONE-1){ printf("ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\nYou can use the 'r' parameter file '%s' which uses option mle=0 to get other results. ",MAXRESULTLINESPONE-1,nresult,rfileres); fprintf(ficlog,"ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\nYou can use the 'r' parameter file '%s' which uses option mle=0 to get other results. ",MAXRESULTLINESPONE-1,nresult,rfileres); @@ -14162,12 +14408,14 @@ Please run with mle=-1 to get a correct /* Tvresult[nres][j] Name of the variable at position j in this resultline */ /* Tresult[nres][j] Value of this variable at position j could be a float if quantitative */ /* We give up with the combinations!! */ - printf("\n j=%d In computing T_ Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=%d cptcovs=%d, cptcoveff=%d Fixed[modelresult[nres][j]]=%d\n", j, nres, j, Dummy[modelresult[nres][j]],nres,j,modelresult[nres][j],cptcovs, cptcoveff,Fixed[modelresult[nres][j]]); /* end if dummy or quanti */ + /* if(debugILK) */ + /* printf("\n j=%d In computing T_ Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=%d cptcovs=%d, cptcoveff=%d Fixed[modelresult[nres][j]]=%d\n", j, nres, j, Dummy[modelresult[nres][j]],nres,j,modelresult[nres][j],cptcovs, cptcoveff,Fixed[modelresult[nres][j]]); /\* end if dummy or quanti *\/ */ if(Dummy[modelresult[nres][j]]==0){/* Dummy variable of the variable in position modelresult in the model corresponding to j in resultline */ - printf("V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ - fprintf(ficlog,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ - fprintf(ficrest,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + /* printf("V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /\* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline *\/ */ /* TinvDoQresult[nres][Name of the variable] */ + printf("V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); /* Output of each value for the combination TKresult[nres], ordered by the covariate values in the resultline */ + fprintf(ficlog,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + fprintf(ficrest,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ if(Fixed[modelresult[nres][j]]==0){ /* Fixed */ printf("fixed ");fprintf(ficlog,"fixed ");fprintf(ficrest,"fixed "); }else{