Summary: Start working on projected prevalences

diff --git a/src/imach.c b/src/imach.c
index 0445b3cfc5db7c8c694a2fe4cd54c4e212d83483..d4647209c1dcb204cdc48d1399cb647c2bdfcec6 100644 (file)
--- a/src/imach.c
+++ b/src/imach.c
@@ -1,6 +1,16 @@
 /* $Id$
   $State$
   $Log$
+  Revision 1.209  2015/11/17 22:12:03  brouard
+  Summary: Adding ftolpl parameter
+  Author: N Brouard
+
+  We had difficulties to get smoothed confidence intervals. It was due
+  to the period prevalence which wasn't computed accurately. The inner
+  parameter ftolpl is now an outer parameter of the .imach parameter
+  file after estepm. If ftolpl is small 1.e-4 and estepm too,
+  computation are long.
+
   Revision 1.208  2015/11/17 14:31:57  brouard
   Summary: temporary
 
@@ -3211,7 +3221,7 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
       k2cpt=0;
       for (i=1; i<=imx; i++) {
        bool=1;
-       if  (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
+       if (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
          for (z1=1; z1<=cptcoveff; z1++)       
             if (covar[Tvaraff[z1]][i]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]){
                 /* Tests if the value of each of the covariates of i is equal to filter j1 */
@@ -3221,7 +3231,7 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
                 j1,z1,nbcode[Tvaraff[z1]][codtabm(j1,z1)],j1);*/
               /* For j1=7 in V1+V2+V3+V4 = 0 1 1 0 and codtabm(7,3)=1 and nbcde[3][?]=1*/
             } 
-       }
+       } /* cptcovn > 0 */
  
        if (bool==1){
          for(m=firstpass; m<=lastpass; m++){
@@ -3235,14 +3245,15 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
                freq[s[m][i]][s[m+1][i]][iagemax+3] += weight[i];
              }
              
-             if ((agev[m][i]>1) && (agev[m][i]< (iagemax+3))) {
+             if ((agev[m][i]>1) && (agev[m][i]< (iagemax+3)) && (anint[m][i]!=9999) && (mint[m][i]!=99)) {
                dateintsum=dateintsum+k2;
                k2cpt++;
+               /* printf("i=%ld dateintmean = %lf dateintsum=%lf k2cpt=%lf k2=%lf\n",i, dateintsum/k2cpt, dateintsum,k2cpt, k2); */
              }
              /*}*/
-         }
-       }
-      } /* end i */
+         } /* end m */
+       } /* end bool */
+      } /* end i = 1 to imx */
        
       /*      fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
       pstamp(ficresp);
@@ -3328,9 +3339,9 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
        if(first==1)
          printf("Others in log...\n");
        fprintf(ficlog,"\n");
-      }
+      } /* end loop i */
       /*}*/
-  }
+  } /* end j1 */
   dateintmean=dateintsum/k2cpt; 
  
   fclose(ficresp);
@@ -4963,6 +4974,11 @@ void printinggnuplot(char fileresu[], char optionfilefiname[], double ageminpar,
     /*#endif */
   m=pow(2,cptcoveff);
 
+  /* Projected Prevalences */
+/* plot "NAGI0w_V1V2_monthlyb2b-proj/F_NAGI0w_V1V2_monthlyb2b-proj.txt" u 6:((($1 == 1) && ($2==0) && ($3==2) &&($4==0))? $7/(1-$13):1/0) t 'p11' w line */
+/* replot ""  u 6:((($1 == 1) && ($2==0) && ($3==2) &&($4==0))? $8/(1-$14):1/0) t 'p21' w line */
+/* replot ""  u 6:((($1 == 1) && ($2==0) && ($3==2) &&($4==0)&&($9!=0))? $9/(1-$15):1/0) t 'p.1' w line */
+
   /* Contribution to likelihood */
   /* Plot the probability implied in the likelihood */
     fprintf(ficgp,"\n# Contributions to the Likelihood, mle >=1. For mle=4 no interpolation, pure matrix products.\n#\n");