]> henry.ined.fr Git - .git/commitdiff
Summary: improvements...
authorN. Brouard <brouard@ined.fr>
Wed, 31 Aug 2022 08:23:16 +0000 (08:23 +0000)
committerN. Brouard <brouard@ined.fr>
Wed, 31 Aug 2022 08:23:16 +0000 (08:23 +0000)
src/imach.c

index 60eb1e939a94cb7f38df2d0e095dd461d1cd80df..d4e416cb9f36f11f8ae77f885296b86ff67724c6 100644 (file)
@@ -1,6 +1,9 @@
 /* $Id$
   $State$
   $Log$
+  Revision 1.334  2022/08/25 09:08:41  brouard
+  Summary: In progress for quantitative
+
   Revision 1.333  2022/08/21 09:10:30  brouard
   * src/imach.c (Module): Version 0.99r33 A lot of changes in
   reassigning covariates: my first idea was that people will always
@@ -1290,13 +1293,13 @@ int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of war
 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 */
-int cptcovn=0; /**< cptcovn decodemodel: number of covariates k of the models excluding age*products =6 and age*age */
+int cptcovn=0; /**< cptcovn decodemodel: number of covariates k of the models excluding age*products =6 and age*age but including products */
 int cptcovt=0; /**< cptcovt: total number of covariates of the model (2) nbocc(+)+1 = 8 excepting constant and age and age*age */
-int cptcovs=0; /**< cptcovs number of simple covariates in the model V2+V1 =2 (dummy or quantit or time varying) */
-int cptcovsnq=0; /**< cptcovsnq number of simple covariates in the model but non quantitative V2+V1 =2 */
+int cptcovs=0; /**< cptcovs number of SIMPLE covariates in the model V2+V1 =2 (dummy or quantit or time varying) */
+int cptcovsnq=0; /**< cptcovsnq number of SIMPLE covariates in the model but non quantitative V2+V1 =2 */
 int cptcovage=0; /**< Number of covariates with age: V3*age only =1 */
 int cptcovprodnoage=0; /**< Number of covariate products without age */   
-int cptcoveff=0; /* Total number of single dummy covariates to vary for printing results (2**cptcoveff combinations of dummies)(computed in tricode as cptcov) */
+int cptcoveff=0; /* Total number of single dummy covariates (fixed or time varying) to vary for printing results (2**cptcoveff combinations of dummies)(computed in tricode as cptcov) */
 int ncovf=0; /* Total number of effective fixed covariates (dummy or quantitative) in the model */
 int ncovv=0; /* Total number of effective (wave) varying covariates (dummy or quantitative) in the model */
 int ncova=0; /* Total number of effective (wave and stepm) varying with age covariates (dummy of quantitative) in the model */
@@ -4169,7 +4172,7 @@ double func( double *x)
 double funcone( double *x)
 {
   /* Same as func but slower because of a lot of printf and if */
-  int i, ii, j, k, mi, d, kk;
+  int i, ii, j, k, mi, d, kk, kf=0;
   int ioffset=0;
   double l, ll[NLSTATEMAX+1], cov[NCOVMAX+1];
   double **out;
@@ -4197,8 +4200,8 @@ double funcone( double *x)
     /* Fixed */
     /* for (k=1; k<=cptcovn;k++) cov[2+nagesqr+k]=covar[Tvar[k]][i]; */
     /* for (k=1; k<=ncoveff;k++){ /\* Simple and product fixed Dummy covariates without age* products *\/ */
-    for (k=1; k<=ncovf;k++){ /* Simple and product fixed covariates without age* products *//* Missing values are set to -1 but should be dropped */
-      cov[ioffset+TvarFind[k]]=covar[Tvar[TvarFind[k]]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V1 is fixed (k=6)*/
+    for (kf=1; kf<=ncovf;kf++){ /* Simple and product fixed covariates without age* products *//* Missing values are set to -1 but should be dropped */
+      cov[ioffset+TvarFind[kf]]=covar[Tvar[TvarFind[kf]]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V1 is fixed (k=6)*/
 /*    cov[ioffset+TvarFind[1]]=covar[Tvar[TvarFind[1]]][i];  */
 /*    cov[2+6]=covar[Tvar[6]][i];  */
 /*    cov[2+6]=covar[2][i]; V2  */
@@ -4298,27 +4301,33 @@ double funcone( 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("Funcone 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],(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])); */
        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);
+       /* printf(" %10.6f\n", -llt); */
       }
-       } /* end of wave */
-} /* end of individual */
-for(k=1,l=0.; k<=nlstate; k++) l += ll[k];
+    } /* end of wave */
+  } /* end of individual */
+  for(k=1,l=0.; k<=nlstate; k++) l += ll[k];
 /* printf("l1=%f l2=%f ",ll[1],ll[2]); */
-l= l*ipmx/sw; /* To get the same order of magnitude as if weight=1 for every body */
-if(globpr==0){ /* First time we count the contributions and weights */
-       gipmx=ipmx;
-       gsw=sw;
-}
+  l= l*ipmx/sw; /* To get the same order of magnitude as if weight=1 for every body */
+  if(globpr==0){ /* First time we count the contributions and weights */
+    gipmx=ipmx;
+    gsw=sw;
+  }
 return -l;
 }
 
@@ -4977,7 +4986,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
   j1=0;
   
   /* j=ncoveff;  /\* Only fixed dummy covariates *\/ */
-  j=cptcoveff;  /* Only dummy covariates used in the model */
+  j=cptcoveff;  /* Only simple dummy covariates used in the model */
   /* j=cptcovn;  /\* Only dummy covariates of the model *\/ */
   if (cptcovn<1) {j=1;ncodemax[1]=1;}
   
@@ -4998,7 +5007,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
 
   /* 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**cptcovn) covariate combination
+   *   Loop on j1(1 to 2**cptcoveff) covariate combination
    *     freq[s1][s2][iage] =0.
    *     Loop on iind
    *       ++freq[s1][s2][iage] weighted
@@ -5023,7 +5032,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
     if(nj==1)
       j=0;  /* First pass for the constant */
     else{
-      j=cptcovs; /* Other passes for the covariate values */
+      j=cptcoveff; /* Other passes for the covariate values number of simple covariates in the model V2+V1 =2 (simple dummy fixed or time varying) */
     }
     first=1;
     for (j1 = 1; j1 <= (int) pow(2,j); j1++){ /* Loop on all dummy covariates combination of the model, ie excluding quantitatives, V4=0, V3=0 for example, fixed or varying covariates */
@@ -5059,13 +5068,16 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
        bool=1;
        if(j !=0){
          if(anyvaryingduminmodel==0){ /* If All fixed covariates */
-           if (cptcovn >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
-             for (z1=1; z1<=cptcovn; z1++) { /* loops on covariates in the model */
+           if (cptcoveff >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
+             for (z1=1; z1<=cptcoveff; z1++) { /* loops on covariates in the model */
                /* if(Tvaraff[z1] ==-20){ */
                /*       /\* sumnew+=cotvar[mw[mi][iind]][z1][iind]; *\/ */
                /* }else  if(Tvaraff[z1] ==-10){ */
                /*       /\* sumnew+=coqvar[z1][iind]; *\/ */
                /* }else  */ /* TODO TODO codtabm(j1,z1) or codtabm(j1,Tvaraff[z1]]z1)*/
+               /* if( iind >=imx-3) printf("Searching error iind=%d Tvaraff[z1]=%d covar[Tvaraff[z1]][iind]=%.f TnsdVar[Tvaraff[z1]]=%d, cptcoveff=%d, cptcovs=%d \n",iind, Tvaraff[z1], covar[Tvaraff[z1]][iind],TnsdVar[Tvaraff[z1]],cptcoveff, cptcovs); */
+               if(Tvaraff[z1]<1 || Tvaraff[z1]>=NCOVMAX)
+                 printf("Error Tvaraff[z1]=%d<1 or >=%d, cptcoveff=%d model=%s\n",Tvaraff[z1],NCOVMAX, cptcoveff, model);
                if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]){ /* for combination j1 of covariates */
                  /* 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 */
@@ -5075,7 +5087,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
                  /* For j1=7 in V1+V2+V3+V4 = 0 1 1 0 and codtabm(7,3)=1 and nbcde[3][?]=1*/
                } /* Onlyf fixed */
              } /* end z1 */
-           } /* cptcovn > 0 */
+           } /* cptcoveff > 0 */
          } /* end any */
        }/* end j==0 */
        if (bool==1){ /* We selected an individual iind satisfying combination j1 (V4=1 V3=0) or all fixed covariates */
@@ -5084,7 +5096,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
            m=mw[mi][iind];
            if(j!=0){
              if(anyvaryingduminmodel==1){ /* Some are varying covariates */
-               for (z1=1; z1<=cptcovn; z1++) {
+               for (z1=1; z1<=cptcoveff; z1++) {
                  if( Fixed[Tmodelind[z1]]==1){
                    iv= Tvar[Tmodelind[z1]]-ncovcol-nqv;
                    if (cotvar[m][iv][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]) /* iv=1 to ntv, right modality. If covariate's 
@@ -5158,9 +5170,9 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
       
       
       /*      fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
-      if(cptcovn==0 && nj==1) /* no covariate and first pass */
+      if(cptcoveff==0 && nj==1) /* no covariate and first pass */
         pstamp(ficresp);
-      if  (cptcovn>0 && j!=0){
+      if  (cptcoveff>0 && j!=0){
         pstamp(ficresp);
        printf( "\n#********** Variable "); 
        fprintf(ficresp, "\n#********** Variable "); 
@@ -5213,14 +5225,17 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
       /* } */
 
       fprintf(ficresphtm,"<table style=\"text-align:center; border: 1px solid\">");
-      if((cptcovn==0 && nj==1)|| nj==2 ) /* no covariate and first pass */
+      if((cptcoveff==0 && nj==1)|| nj==2 ) /* no covariate and first pass */
         fprintf(ficresp, " Age");
-      if(nj==2) for (z1=1; z1<=cptcovn; z1++) fprintf(ficresp, " V%d=%d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]);
+      if(nj==2) for (z1=1; z1<=cptcoveff; z1++) {
+         printf(" V%d=%d, z1=%d, Tvaraff[z1]=%d, j1=%d, TnsdVar[Tvaraff[%d]]=%d |",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])], z1, Tvaraff[z1], j1,z1,TnsdVar[Tvaraff[z1]]);
+         fprintf(ficresp, " V%d=%d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]);
+       }
       for(i=1; i<=nlstate;i++) {
-       if((cptcovn==0 && nj==1)|| nj==2 ) fprintf(ficresp," 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, "<th>Age</th><th>Prev(%d)</th><th>N(%d)</th><th>N</th>",i,i);
       }
-      if((cptcovn==0 && nj==1)|| nj==2 ) fprintf(ficresp, "\n");
+      if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp, "\n");
       fprintf(ficresphtm, "\n");
       
       /* Header of frequency table by age */
@@ -5288,14 +5303,14 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
        }
        
        /* Writing ficresp */
-       if(cptcovn==0 && nj==1){ /* no covariate and first pass */
+       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<=cptcovn; z1++) fprintf(ficresp, " %d %d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]);
+            for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresp, " %d %d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]);
           }
        }
        for(s1=1; s1 <=nlstate ; s1++){
@@ -5310,7 +5325,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
          }
          if( iage <= iagemax){
            if(pos>=1.e-5){
-             if(cptcovn==0 && nj==1){ /* no covariate and first pass */
+             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);
@@ -5319,7 +5334,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
              /*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((cptcovn==0 && nj==1)|| nj==2 ) fprintf(ficresp," NaNq %.0f %.0f",prop[s1][iage],pospropta);
+             if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp," NaNq %.0f %.0f",prop[s1][iage],pospropta);
              fprintf(ficresphtm,"<th>%d</th><td>NaNq</td><td>%.0f</td><td>%.0f</td>",iage, prop[s1][iage],pospropta);
            }
          }
@@ -5345,7 +5360,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
        }
        fprintf(ficresphtmfr,"</tr>\n ");
        fprintf(ficresphtm,"</tr>\n");
-       if((cptcovn==0 && nj==1)|| nj==2 ) {
+       if((cptcoveff==0 && nj==1)|| nj==2 ) {
          if(iage <= iagemax)
            fprintf(ficresp,"\n");
        }
@@ -5603,7 +5618,7 @@ void prevalence(double ***probs, double agemin, double agemax, int **s, double *
   if (cptcovn<1) {j=1;ncodemax[1]=1;}
   
   first=0;
-  for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of covariate */
+  for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of simple dummy covariates */
     for (i=1; i<=nlstate; i++)  
       for(iage=iagemin-AGEMARGE; iage <= iagemax+4+AGEMARGE; iage++)
        prop[i][iage]=0.0;
@@ -5927,7 +5942,7 @@ void  concatwav(int wav[], int **dh, int **bh,  int **mw, int **s, double *agedc
        nbcode[k][j]=0; /* Valgrind */
 
    /* Loop on covariates without age and products and no quantitative variable */
-   for (k=1; k<=cptcovt; k++) { /* From model V1 + V2*age + V3 + V3*V4 keeps V1 + V3 = 2 only */
+   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;
      if(Dummy[k]==0 && Typevar[k] !=1){ /* Dummy covariate and not age product */ 
        switch(Fixed[k]) {
@@ -6027,6 +6042,10 @@ void  concatwav(int wav[], int **dh, int **bh,  int **mw, int **s, double *agedc
      } /* end dummy test */
      if(Dummy[k]==1 && Typevar[k] !=1){ /* 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);
+          exit(1);
+        }
         if(isnan(covar[Tvar[k]][i])){
           printf("ERROR, IMaCh doesn't treat fixed quantitative covariate with missing values V%d=., individual %d will be skipped.\n",Tvar[k],i);
           fprintf(ficlog,"ERROR, currently IMaCh doesn't treat covariate with missing values V%d=., individual %d will be skipped.\n",Tvar[k],i);
@@ -6034,7 +6053,7 @@ void  concatwav(int wav[], int **dh, int **bh,  int **mw, int **s, double *agedc
           exit(1);
          }
        }
-     }
+     } /* end Quanti */
    } /* end of loop on model-covariate k. nbcode[Tvark][1]=-1, nbcode[Tvark][1]=0 and nbcode[Tvark][2]=1 sets the value of covariate k*/  
   
    for (k=-1; k< maxncov; k++) Ndum[k]=0; 
@@ -6048,13 +6067,22 @@ void  concatwav(int wav[], int **dh, int **bh,  int **mw, int **s, double *agedc
   
    ij=0;
    /* for (i=0; i<=  maxncov-1; i++) { /\* modmaxcovj is unknown here. Only Ndum[2(V2),3(age*V3), 5(V3*V2) 6(V1*V4) *\/ */
-   for (k=1; k<=  cptcovt; k++) { /* modmaxcovj is unknown here. Only Ndum[2(V2),3(age*V3), 5(V3*V2) 6(V1*V4) */
+   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 */
+     /* modmaxcovj is unknown here. Only Ndum[2(V2),3(age*V3), 5(V3*V2) 6(V1*V4) */
      /*printf("Ndum[%d]=%d\n",i, Ndum[i]);*/
      /* if((Ndum[i]!=0) && (i<=ncovcol)){  /\* Tvar[i] <= ncovmodel ? *\/ */
-     if(Ndum[Tvar[k]]!=0 && Dummy[k] == 0 && Typevar[k]==0){  /* Only Dummy and non empty in the model */
+     if(Ndum[Tvar[k]]!=0 && Dummy[k] == 0 && Typevar[k]==0){  /* Only Dummy simple and non empty in the model */
+       /* Typevar[k] =0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product */
+       /* Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product*/
        /* If product not in single variable we don't print results */
        /*printf("diff Ndum[%d]=%d\n",i, Ndum[i]);*/
-       ++ij;/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, */
+       ++ij;/*    V5 + V4 + V3 + V4*V3 + V5*age + V2 +  V1*V2 + V1*age + V1, *//* V5 quanti, V2 quanti, V4, V3, V1 dummies */
+       /* k=       1    2   3     4       5       6      7       8        9  */
+       /* Tvar[k]= 5    4    3    6       5       2      7       1        1  */
+       /* ij            1    2                                            3  */  
+       /* Tvaraff[ij]=  4    3                                            1  */
+       /* Tmodelind[ij]=2    3                                            9  */
+       /* TmodelInvind[ij]=2 1                                            1  */
        Tvaraff[ij]=Tvar[k]; /* For printing combination *//* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, Tvar {5, 4, 3, 6, 5, 2, 7, 1, 1} Tvaraff={4, 3, 1} V4, V3, V1*/
        Tmodelind[ij]=k; /* Tmodelind: index in model of dummies Tmodelind[1]=2 V4: pos=2; V3: pos=3, V1=9 {2, 3, 9, ?, ?,} */
        TmodelInvind[ij]=Tvar[k]- ncovcol-nqv; /* Inverse TmodelInvind[2=V4]=2 second dummy varying cov (V4)4-1-1 {0, 2, 1, } TmodelInvind[3]=1 */
@@ -6070,7 +6098,7 @@ void  concatwav(int wav[], int **dh, int **bh,  int **mw, int **s, double *agedc
    } /* Tvaraff[1]@5 {3, 4, -20, 0, 0} Very strange */
    /* ij--; */
    /* cptcoveff=ij; /\*Number of total covariates*\/ */
-   *cptcov=ij; /* cptcov= Number of total real effective covariates: effective (used as cptcoveff in other functions)
+   *cptcov=ij; /* cptcov= Number of total real effective simple dummies (fixed or time  arying) effective (used as cptcoveff in other functions)
                * because they can be excluded from the model and real
                * if in the model but excluded because missing values, but how to get k from ij?*/
    for(j=ij+1; j<= cptcovt; j++){
@@ -10704,9 +10732,17 @@ int decodemodel( char model[], int lastobs)
            Tvar[k]=ncovcol+nqv+ntv+nqtv+k1; /* For model-covariate k tells which data-covariate to use but
                                                because this model-covariate is a construction we invent a new column
                                                which is after existing variables ncovcol+nqv+ntv+nqtv + k1
-                                               If already ncovcol=4 and model=V2 + V1 +V1*V4 +age*V3 +V3*V2
+                                               If already ncovcol=4 and model= V2 + V1 + V1*V4 + age*V3 + V3*V2
                                                thus after V4 we invent V5 and V6 because age*V3 will be computed in 4
                                                Tvar[3=V1*V4]=4+1=5 Tvar[5=V3*V2]=4 + 2= 6, Tvar[4=age*V3]=4 etc */
+           /* Please remark that the new variables are model dependent */
+           /* If we have 4 variable but the model uses only 3, like in
+            * model= V1 + age*V1 + V2 + V3 + age*V2 + age*V3 + V1*V2 + V1*V3
+            *  k=     1     2       3   4     5        6        7       8
+            * Tvar[k]=1     1       2   3     2        3       (5       6) (and not 4 5 because of V4 missing)
+            * Tage[kk]    [1]= 2           [2]=5      [3]=6                  kk=1 to cptcovage=3
+            * Tvar[Tage[kk]][1]=2          [2]=2      [3]=3
+            */
            Typevar[k]=2;  /* 2 for double fixed dummy covariates */
            cutl(strc,strb,strd,'V'); /* strd was Vm, strc is m */
            Tprod[k1]=k;  /* Tprod[1]=3(=V1*V4) for V2+V1+V1*V4+age*V3+V3*V2  */
@@ -12428,6 +12464,7 @@ Please run with mle=-1 to get a correct covariance matrix.\n",optionfile,numline
   Tvar=ivector(1,NCOVMAX); /* Was 15 changed to NCOVMAX. */
   TvarsDind=ivector(1,NCOVMAX); /*  */
   TnsdVar=ivector(1,NCOVMAX); /*  */
+    /* for(i=1; i<=NCOVMAX;i++) TnsdVar[i]=3; */
   TvarsD=ivector(1,NCOVMAX); /*  */
   TvarsQind=ivector(1,NCOVMAX); /*  */
   TvarsQ=ivector(1,NCOVMAX); /*  */
@@ -12550,7 +12587,7 @@ Please run with mle=-1 to get a correct covariance matrix.\n",optionfile,numline
   Ndum =ivector(-1,NCOVMAX);  
   cptcoveff=0;
   if (ncovmodel-nagesqr > 2 ){ /* That is if covariate other than cst, age and age*age */
-    tricode(&cptcoveff,Tvar,nbcode,imx, Ndum); /**< Fills nbcode[Tvar[j]][l]; */
+    tricode(&cptcoveff,Tvar,nbcode,imx, Ndum); /**< Fills nbcode[Tvar[j]][l]; as well as calculate cptcoveff or number of total effective dummy covariates*/
   }
   
   ncovcombmax=pow(2,cptcoveff);
@@ -12689,11 +12726,18 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
          optionfilehtmcov,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model);
   }
 
-  fprintf(fichtm,"<html><head>\n<head>\n<meta charset=\"utf-8\"/><meta http-equiv=\"Content-Type\" content=\"text/html; charset=utf-8\" />\n<title>IMaCh %s</title></head>\n <body><font size=\"7\"><a href=http:/euroreves.ined.fr/imach>IMaCh for Interpolated Markov Chain</a> </font><br>\n<font size=\"3\">Sponsored by Copyright (C)  2002-2015 <a href=http://www.ined.fr>INED</a>-EUROREVES-Institut de longévité-2013-2022-Japan Society for the Promotion of Sciences 日本学術振興会 (<a href=https://www.jsps.go.jp/english/e-grants/>Grant-in-Aid for Scientific Research 25293121</a>) - <a href=https://software.intel.com/en-us>Intel Software 2015-2018</a></font><br>  \
-<hr size=\"2\" color=\"#EC5E5E\"> \n\
+  fprintf(fichtm,"<html><head>\n<meta charset=\"utf-8\"/><meta http-equiv=\"Content-Type\" content=\"text/html; charset=utf-8\" />\n\
+<title>IMaCh %s</title></head>\n\
+ <body><font size=\"7\"><a href=http:/euroreves.ined.fr/imach>IMaCh for Interpolated Markov Chain</a> </font><br>\n\
+<font size=\"3\">Sponsored by Copyright (C)  2002-2015 <a href=http://www.ined.fr>INED</a>\
+-EUROREVES-Institut de longévité-2013-2022-Japan Society for the Promotion of Sciences 日本学術振興会 \
+(<a href=https://www.jsps.go.jp/english/e-grants/>Grant-in-Aid for Scientific Research 25293121</a>) - \
+<a href=https://software.intel.com/en-us>Intel Software 2015-2018</a></font><br> \n", optionfilehtm);
+  
+  fprintf(fichtm,"<hr size=\"2\" color=\"#EC5E5E\"> \n\
 <font size=\"2\">IMaCh-%s <br> %s</font> \
 <hr size=\"2\" color=\"#EC5E5E\"> \n\
-Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s<br>\n\
+This file: <a href=\"%s\">%s</a>Title=%s <br>Datafile=<a href=\"%s\">%s</a> Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s<br>\n\
 \n\
 <hr  size=\"2\" color=\"#EC5E5E\">\
  <ul><li><h4>Parameter files</h4>\n\
@@ -12702,7 +12746,7 @@ Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age
  - Log file of the run: <a href=\"%s\">%s</a><br>\n\
  - Gnuplot file name: <a href=\"%s\">%s</a><br>\n\
  - Date and time at start: %s</ul>\n",\
-         optionfilehtm,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model,\
+         version,fullversion,optionfilehtm,optionfilehtm,title,datafile,datafile,firstpass,lastpass,stepm, weightopt, model, \
          optionfilefiname,optionfilext,optionfilefiname,optionfilext,\
          fileres,fileres,\
          filelog,filelog,optionfilegnuplot,optionfilegnuplot,strstart);
@@ -13016,6 +13060,7 @@ Please run with mle=-1 to get a correct covariance matrix.\n",ageminpar,agemaxpa
     globpr=1; /* again, to print the individual contributions using computed gpimx and gsw */
     likelione(ficres, p, npar, nlstate, &globpr, &ipmx, &sw, &fretone, funcone); /* Prints the contributions to the likelihood */
     printf("Second Likeli=%12.6f ipmx=%ld sw=%12.6f",fretone,ipmx,sw);
+          /* exit(0); */
     for (k=1; k<=npar;k++)
       printf(" %d %8.5f",k,p[k]);
     printf("\n");
@@ -13810,7 +13855,7 @@ Please run with mle=-1 to get a correct covariance matrix.\n",ageminpar,agemaxpa
          fprintf(ficlog,"Error in computing T_ Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=%d cptcovs=%d, cptcoveff=%d \n", nres, j, Dummy[modelresult[nres][j]],nres,j,modelresult[nres][j],cptcovs, cptcoveff);  /* end if dummy  or quanti */
          exit(1);
        }
-      }
+      } /* End loop for each variable in the resultline */
       /* for (j=1; j<= nsq; j++){ /\* For each selected (single) quantitative value *\/ */
       /*       printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); /\* Wrong j is not in the equation model *\/ */
       /*       fprintf(ficrest," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */
@@ -13837,7 +13882,8 @@ Please run with mle=-1 to get a correct covariance matrix.\n",ageminpar,agemaxpa
       fprintf(ficresvij,"\n#****** ");
       /* pstamp(ficresvij); */
       for(j=1;j<=cptcoveff;j++) 
-       fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[TnsdVar[Tvaraff[j]]])]);
+       fprintf(ficresvij,"V%d=%d ",Tvresult[nres][j],Tresult[nres][j]);
+       /* fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[TnsdVar[Tvaraff[j]]])]); */
       for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
        /* fprintf(ficresvij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]); /\* To solve *\/ */
        fprintf(ficresvij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); /* Solved */
@@ -13871,7 +13917,7 @@ Please run with mle=-1 to get a correct covariance matrix.\n",ageminpar,agemaxpa
          fprintf(ficrest,"the age specific prevalence observed (cross-sectionally) in the population i.e cross-sectionally\n in each health state (popbased=1) (mobilav=%d)\n",mobilav);
        else
          fprintf(ficrest,"the age specific forward period (stable) prevalences in each health state \n");
-       fprintf(ficrest,"# Age popbased mobilav e.. (std) ");
+       fprintf(ficrest,"# Age popbased mobilav e.. (std) "); /* Adding covariate values? */
        for (i=1;i<=nlstate;i++) fprintf(ficrest,"e.%d (std) ",i);
        fprintf(ficrest,"\n");
        /* printf("Which p?\n"); for(i=1;i<=npar;i++)printf("p[i=%d]=%lf,",i,p[i]);printf("\n"); */
@@ -13918,12 +13964,12 @@ Please run with mle=-1 to get a correct covariance matrix.\n",ageminpar,agemaxpa
       printf("done selection\n");fflush(stdout);
       fprintf(ficlog,"done selection\n");fflush(ficlog);
       
-    } /* End k selection */
+    } /* End k selection or end covariate selection for nres */
 
     printf("done State-specific expectancies\n");fflush(stdout);
     fprintf(ficlog,"done State-specific expectancies\n");fflush(ficlog);
 
-    /* variance-covariance of forward period prevalence*/
+    /* variance-covariance of forward period prevalence */
     varprlim(fileresu, nresult, mobaverage, mobilavproj, bage, fage, prlim, &ncvyear, ftolpl, p, matcov, delti, stepm, cptcoveff);