--- imach/src/imach.c 2016/08/26 09:20:19 1.237
+++ imach/src/imach.c 2016/08/29 17:17:25 1.241
@@ -1,6 +1,20 @@
-/* $Id: imach.c,v 1.237 2016/08/26 09:20:19 brouard Exp $
+/* $Id: imach.c,v 1.241 2016/08/29 17:17:25 brouard Exp $
$State: Exp $
$Log: imach.c,v $
+ Revision 1.241 2016/08/29 17:17:25 brouard
+ Summary: gnuplot problem in Back projection to fix
+
+ Revision 1.240 2016/08/29 07:53:18 brouard
+ Summary: Better
+
+ Revision 1.239 2016/08/26 15:51:03 brouard
+ Summary: Improvement in Powell output in order to copy and paste
+
+ Author:
+
+ Revision 1.238 2016/08/26 14:23:35 brouard
+ Summary: Starting tests of 0.99
+
Revision 1.237 2016/08/26 09:20:19 brouard
Summary: to valgrind
@@ -911,12 +925,12 @@ typedef struct {
#define ODIRSEPARATOR '\\'
#endif
-/* $Id: imach.c,v 1.237 2016/08/26 09:20:19 brouard Exp $ */
+/* $Id: imach.c,v 1.241 2016/08/29 17:17:25 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.237 $ $Date: 2016/08/26 09:20:19 $";
+char fullversion[]="$Revision: 1.241 $ $Date: 2016/08/29 17:17:25 $";
char strstart[80];
char optionfilext[10], optionfilefiname[FILENAMELENGTH];
int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */
@@ -1029,7 +1043,8 @@ double dval;
#define FTOL 1.0e-10
#define NRANSI
-#define ITMAX 200
+#define ITMAX 200
+#define ITPOWMAX 20 /* This is now multiplied by the number of parameters */
#define TOL 2.0e-4
@@ -1145,6 +1160,8 @@ double *Tvalsel; /**< Selected modality
int *Typevar; /**< 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product */
int *Fixed; /** Fixed[k] 0=fixed, 1 varying, 2 fixed with age product, 3 varying with age product */
int *Dummy; /** Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product */
+int *DummyV; /** Dummy[v] 0=dummy (0 1), 1 quantitative */
+int *FixedV; /** FixedV[v] 0 fixed, 1 varying */
int *Tage;
int anyvaryingduminmodel=0; /**< Any varying dummy in Model=1 yes, 0 no, to avoid a loop on waves in freq */
int *Tmodelind; /** Tmodelind[Tvaraff[3]]=9 for V1 position,Tvaraff[1]@9={4, 3, 1, 0, 0, 0, 0, 0, 0}, model=V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1*/
@@ -1154,11 +1171,10 @@ int *Ndum; /** Freq of modality (tricode
/* int **codtab;*/ /**< codtab=imatrix(1,100,1,10); */
int **Tvard;
int *Tprod;/**< Gives the k position of the k1 product */
+/* Tprod[k1=1]=3(=V1*V4) for V2+V1+V1*V4+age*V3 */
int *Tposprod; /**< Gives the k1 product from the k position */
-/* Tprod[k1=1]=3(=V1*V4) for V2+V1+V1*V4+age*V3
- if V2+V1+V1*V4+age*V3+V3*V2 TProd[k1=2]=5 (V3*V2)
- Tposprod[k]=k1 , Tposprod[3]=1, Tposprod[5]=2
-*/
+ /* if V2+V1+V1*V4+age*V3+V3*V2 TProd[k1=2]=5 (V3*V2) */
+ /* Tposprod[k]=k1 , Tposprod[3]=1, Tposprod[5(V3*V2)]=2 (2nd product without age) */
int cptcovprod, *Tvaraff, *invalidvarcomb;
double *lsurv, *lpop, *tpop;
@@ -2057,10 +2073,10 @@ void powell(double p[], double **xi, int
void linmin(double p[], double xi[], int n, double *fret,
double (*func)(double []));
#else
- void linmin(double p[], double xi[], int n, double *fret,
- double (*func)(double []),int *flat);
+ void linmin(double p[], double xi[], int n, double *fret,
+ double (*func)(double []),int *flat);
#endif
- int i,ibig,j;
+ int i,ibig,j,jk,k;
double del,t,*pt,*ptt,*xit;
double directest;
double fp,fptt;
@@ -2092,31 +2108,64 @@ void powell(double p[], double **xi, int
fprintf(ficlog,"\nPowell iter=%d -2*LL=%.12f %ld sec. %ld sec.",*iter,*fret,rcurr_time-rlast_time, rcurr_time-rstart_time); fflush(ficlog);
/* fprintf(ficrespow,"%d %.12f %ld",*iter,*fret,curr_time.tm_sec-start_time.tm_sec); */
for (i=1;i<=n;i++) {
- printf(" %d %.12f",i, p[i]);
- fprintf(ficlog," %d %.12lf",i, p[i]);
fprintf(ficrespow," %.12lf", p[i]);
}
+ fprintf(ficrespow,"\n");fflush(ficrespow);
+ printf("\n#model= 1 + age ");
+ fprintf(ficlog,"\n#model= 1 + age ");
+ if(nagesqr==1){
+ printf(" + age*age ");
+ fprintf(ficlog," + age*age ");
+ }
+ for(j=1;j <=ncovmodel-2;j++){
+ if(Typevar[j]==0) {
+ printf(" + V%d ",Tvar[j]);
+ fprintf(ficlog," + V%d ",Tvar[j]);
+ }else if(Typevar[j]==1) {
+ printf(" + V%d*age ",Tvar[j]);
+ fprintf(ficlog," + V%d*age ",Tvar[j]);
+ }else if(Typevar[j]==2) {
+ printf(" + V%d*V%d ",Tvard[Tposprod[j]][1],Tvard[Tposprod[j]][2]);
+ fprintf(ficlog," + V%d*V%d ",Tvard[Tposprod[j]][1],Tvard[Tposprod[j]][2]);
+ }
+ }
printf("\n");
+/* printf("12 47.0114589 0.0154322 33.2424412 0.3279905 2.3731903 */
+/* 13 -21.5392400 0.1118147 1.2680506 1.2973408 -1.0663662 */
fprintf(ficlog,"\n");
- fprintf(ficrespow,"\n");fflush(ficrespow);
- if(*iter <=3){
+ for(i=1,jk=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(j=1; j <=ncovmodel; j++){
+ printf("%12.7f ",p[jk]);
+ fprintf(ficlog,"%12.7f ",p[jk]);
+ jk++;
+ }
+ printf("\n");
+ fprintf(ficlog,"\n");
+ }
+ }
+ }
+ if(*iter <=3 && *iter >1){
tml = *localtime(&rcurr_time);
strcpy(strcurr,asctime(&tml));
rforecast_time=rcurr_time;
itmp = strlen(strcurr);
if(strcurr[itmp-1]=='\n') /* Windows outputs with a new line */
- strcurr[itmp-1]='\0';
+ strcurr[itmp-1]='\0';
printf("\nConsidering the time needed for the last iteration #%d: %ld seconds,\n",*iter,rcurr_time-rlast_time);
fprintf(ficlog,"\nConsidering the time needed for this last iteration #%d: %ld seconds,\n",*iter,rcurr_time-rlast_time);
for(niterf=10;niterf<=30;niterf+=10){
- rforecast_time=rcurr_time+(niterf-*iter)*(rcurr_time-rlast_time);
- forecast_time = *localtime(&rforecast_time);
- strcpy(strfor,asctime(&forecast_time));
- itmp = strlen(strfor);
- if(strfor[itmp-1]=='\n')
- strfor[itmp-1]='\0';
- printf(" - if your program needs %d iterations to converge, convergence will be \n reached in %s i.e.\n on %s (current time is %s);\n",niterf, asc_diff_time(rforecast_time-rcurr_time,tmpout),strfor,strcurr);
- fprintf(ficlog," - if your program needs %d iterations to converge, convergence will be \n reached in %s i.e.\n on %s (current time is %s);\n",niterf, asc_diff_time(rforecast_time-rcurr_time,tmpout),strfor,strcurr);
+ rforecast_time=rcurr_time+(niterf-*iter)*(rcurr_time-rlast_time);
+ forecast_time = *localtime(&rforecast_time);
+ strcpy(strfor,asctime(&forecast_time));
+ itmp = strlen(strfor);
+ if(strfor[itmp-1]=='\n')
+ strfor[itmp-1]='\0';
+ printf(" - if your program needs %d iterations to converge, convergence will be \n reached in %s i.e.\n on %s (current time is %s);\n",niterf, asc_diff_time(rforecast_time-rcurr_time,tmpout),strfor,strcurr);
+ fprintf(ficlog," - if your program needs %d iterations to converge, convergence will be \n reached in %s i.e.\n on %s (current time is %s);\n",niterf, asc_diff_time(rforecast_time-rcurr_time,tmpout),strfor,strcurr);
}
}
for (i=1;i<=n;i++) { /* For each direction i */
@@ -2216,7 +2265,7 @@ void powell(double p[], double **xi, int
free_vector(pt,1,n);
return;
} /* enough precision */
- if (*iter == ITMAX) nrerror("powell exceeding maximum iterations.");
+ if (*iter == ITMAX*n) nrerror("powell exceeding maximum iterations.");
for (j=1;j<=n;j++) { /* Computes the extrapolated point P_0 + 2 (P_n-P_0) */
ptt[j]=2.0*p[j]-pt[j];
xit[j]=p[j]-pt[j];
@@ -4090,7 +4139,7 @@ void freqsummary(char fileres[], int ia
fprintf(ficlog,"Problem with prevalence resultfile: %s\n", fileresp);
exit(0);
}
-
+
strcpy(fileresphtm,subdirfext(optionfilefiname,"PHTM_",".htm"));
if((ficresphtm=fopen(fileresphtm,"w"))==NULL) {
printf("Problem with prevalence HTM resultfile '%s' with errno='%s'\n",fileresphtm,strerror(errno));
@@ -4100,55 +4149,54 @@ void freqsummary(char fileres[], int ia
}
else{
fprintf(ficresphtm,"
\nIMaCh PHTM_ %s\n %s
%s \
-
\n\
+
\n \
Title=%s
Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s
\n",\
fileresphtm,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model);
}
fprintf(ficresphtm,"Current page is file %s
\n\nFrequencies and prevalence by age at begin of transition and dummy covariate value at beginning of transition
\n",fileresphtm, fileresphtm);
-
+
strcpy(fileresphtmfr,subdirfext(optionfilefiname,"PHTMFR_",".htm"));
if((ficresphtmfr=fopen(fileresphtmfr,"w"))==NULL) {
printf("Problem with frequency table HTM resultfile '%s' with errno='%s'\n",fileresphtmfr,strerror(errno));
fprintf(ficlog,"Problem with frequency table HTM resultfile '%s' with errno='%s'\n",fileresphtmfr,strerror(errno));
fflush(ficlog);
exit(70);
- }
- else{
+ } else{
fprintf(ficresphtmfr,"\nIMaCh PHTM_Frequency table %s\n %s
%s \
-
\n\
+
\n \
Title=%s
Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s
\n",\
fileresphtmfr,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model);
}
- fprintf(ficresphtmfr,"Current page is file %s
\n\nFrequencies of all effective transitions by age at begin of transition
Unknown status is -1
\n",fileresphtmfr, fileresphtmfr);
-
+ fprintf(ficresphtmfr,"Current page is file %s
\n\nFrequencies of all effective transitions of the model, by age at begin of transition, and covariate value at the begin of transition (if the covariate is a varying covariate)
Unknown status is -1
\n",fileresphtmfr, fileresphtmfr);
+
freq= ma3x(-5,nlstate+ndeath,-5,nlstate+ndeath,iagemin-AGEMARGE,iagemax+3+AGEMARGE);
j1=0;
/* j=ncoveff; /\* Only fixed dummy covariates *\/ */
j=cptcoveff; /* Only dummy covariates of the model */
if (cptcovn<1) {j=1;ncodemax[1]=1;}
-
+
first=1;
-
+
/* Detects if a combination j1 is empty: for a multinomial variable like 3 education levels:
reference=low_education V1=0,V2=0
med_educ V1=1 V2=0,
high_educ V1=0 V2=1
Then V1=1 and V2=1 is a noisy combination that we want to exclude for the list 2**cptcoveff
*/
-
+
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 */
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(m=iagemin; m <= iagemax+3; m++)
- freq[i][jk][m]=0;
-
+ for(m=iagemin; m <= iagemax+3; m++)
+ freq[i][jk][m]=0;
+
for (i=1; i<=nlstate; i++) {
for(m=iagemin; m <= iagemax+3; m++)
- prop[i][m]=0;
+ prop[i][m]=0;
posprop[i]=0;
pospropt[i]=0;
}
@@ -4158,7 +4206,7 @@ Title=%s
Datafile=%s Firstpass=%d La
/* meanqt[m][z1]=0.; */
/* } */
/* } */
-
+
dateintsum=0;
k2cpt=0;
/* For that combination of covariate j1, we count and print the frequencies in one pass */
@@ -4237,35 +4285,41 @@ Title=%s
Datafile=%s Firstpass=%d La
} /* end iind = 1 to imx */
/* prop[s][age] is feeded for any initial and valid live state as well as
freq[s1][s2][age] at single age of beginning the transition, for a combination j1 */
-
-
+
+
/* fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
pstamp(ficresp);
- /* if (ncoveff>0) { */
- if (cptcoveff>0) {
+ if (cptcoveff>0){
fprintf(ficresp, "\n#********** Variable ");
fprintf(ficresphtm, "\n
********** Variable ");
fprintf(ficresphtmfr, "\n
********** Variable ");
+ fprintf(ficlog, "\n#********** Variable ");
for (z1=1; z1<=cptcoveff; z1++){
- fprintf(ficresp, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
- fprintf(ficresphtm, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
- fprintf(ficresphtmfr, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ if(DummyV[z1]){
+ fprintf(ficresp, "V%d (fixed)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresphtm, "V%d (fixed)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresphtmfr, "V%d (fixed)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficlog, "V%d (fixed)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ }else{
+ fprintf(ficresp, "V%d(varying)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresphtm, "V%d(varying)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresphtmfr, "V%d(varying)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficlog, "V%d(varying)=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ }
}
fprintf(ficresp, "**********\n#");
fprintf(ficresphtm, "**********
\n");
fprintf(ficresphtmfr, "**********
\n");
- fprintf(ficlog, "\n#********** Variable ");
- for (z1=1; z1<=cptcoveff; z1++) fprintf(ficlog, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
fprintf(ficlog, "**********\n");
}
fprintf(ficresphtm,"");
for(i=1; i<=nlstate;i++) {
- fprintf(ficresp, " Age Prev(%d) N(%d) N",i,i);
+ fprintf(ficresp, " Age Prev(%d) N(%d) N ",i,i);
fprintf(ficresphtm, "Age | Prev(%d) | N(%d) | N | ",i,i);
}
fprintf(ficresp, "\n");
fprintf(ficresphtm, "\n");
-
+
/* Header of frequency table by age */
fprintf(ficresphtmfr,"");
fprintf(ficresphtmfr,"Age | ");
@@ -4276,115 +4330,115 @@ Title=%s
Datafile=%s Firstpass=%d La
}
}
fprintf(ficresphtmfr, "\n");
-
+
/* For each age */
for(iage=iagemin; iage <= iagemax+3; iage++){
fprintf(ficresphtm,"");
if(iage==iagemax+1){
- fprintf(ficlog,"1");
- fprintf(ficresphtmfr,"
0 | ");
+ fprintf(ficlog,"1");
+ fprintf(ficresphtmfr,"
---|
0 | ");
}else if(iage==iagemax+2){
- fprintf(ficlog,"0");
- fprintf(ficresphtmfr,"
---|
Unknown | ");
+ fprintf(ficlog,"0");
+ fprintf(ficresphtmfr,"
---|
Unknown | ");
}else if(iage==iagemax+3){
- fprintf(ficlog,"Total");
- fprintf(ficresphtmfr,"
---|
Total | ");
+ fprintf(ficlog,"Total");
+ fprintf(ficresphtmfr,"
---|
Total | ");
}else{
- if(first==1){
- first=0;
- printf("See log file for details...\n");
- }
- fprintf(ficresphtmfr,"
---|
%d | ",iage);
- fprintf(ficlog,"Age %d", iage);
+ if(first==1){
+ first=0;
+ printf("See log file for details...\n");
+ }
+ fprintf(ficresphtmfr,"
---|
%d | ",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(m=-1, pp[jk]=0; m <=nlstate+ndeath ; m++)
+ pp[jk] += freq[jk][m][iage];
}
for(jk=1; jk <=nlstate ; jk++){
- for(m=-1, pos=0; m <=0 ; m++)
- pos += freq[jk][m][iage];
- if(pp[jk]>=1.e-10){
- if(first==1){
- printf(" %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
- }
- fprintf(ficlog," %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
- }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);
- }
+ for(m=-1, pos=0; m <=0 ; m++)
+ pos += freq[jk][m][iage];
+ if(pp[jk]>=1.e-10){
+ if(first==1){
+ printf(" %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
+ }
+ fprintf(ficlog," %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
+ }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);
+ }
}
-
+
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];
+ /* 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(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] */
+ 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(jk=1; jk <=nlstate ; jk++){
- 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);
- }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);
- }
- if( iage <= iagemax){
- if(pos>=1.e-5){
- fprintf(ficresp," %d %.5f %.0f %.0f",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta);
- fprintf(ficresphtm,"%d | %.5f | %.0f | %.0f | ",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,"%d | NaNq | %.0f | %.0f | ",iage, prop[jk][iage],pospropta);
- }
- }
- pospropt[jk] +=posprop[jk];
+ 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);
+ }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);
+ }
+ if( iage <= iagemax){
+ if(pos>=1.e-5){
+ fprintf(ficresp," %d %.5f %.0f %.0f",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta);
+ fprintf(ficresphtm,"%d | %.5f | %.0f | %.0f | ",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,"%d | NaNq | %.0f | %.0f | ",iage, prop[jk][iage],pospropta);
+ }
+ }
+ pospropt[jk] +=posprop[jk];
} /* end loop jk */
/* pospropt=0.; */
for(jk=-1; jk <=nlstate+ndeath; jk++){
- for(m=-1; m <=nlstate+ndeath; m++){
- if(freq[jk][m][iage] !=0 ) { /* minimizing output */
- if(first==1){
- printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]);
- }
- fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iage]);
- }
- if(jk!=0 && m!=0)
- fprintf(ficresphtmfr,"%.0f | ",freq[jk][m][iage]);
- }
+ for(m=-1; m <=nlstate+ndeath; m++){
+ if(freq[jk][m][iage] !=0 ) { /* minimizing output */
+ if(first==1){
+ printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]);
+ }
+ fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iage]);
+ }
+ if(jk!=0 && m!=0)
+ fprintf(ficresphtmfr,"%.0f | ",freq[jk][m][iage]);
+ }
} /* end loop jk */
posproptt=0.;
for(jk=1; jk <=nlstate; jk++){
- posproptt += pospropt[jk];
+ posproptt += pospropt[jk];
}
fprintf(ficresphtmfr,"
\n ");
if(iage <= iagemax){
- fprintf(ficresp,"\n");
- fprintf(ficresphtm,"\n");
+ fprintf(ficresp,"\n");
+ fprintf(ficresphtm,"\n");
}
if(first==1)
- printf("Others in log...\n");
+ printf("Others in log...\n");
fprintf(ficlog,"\n");
} /* end loop age iage */
fprintf(ficresphtm,"Tot | ");
for(jk=1; jk <=nlstate ; jk++){
if(posproptt < 1.e-5){
- fprintf(ficresphtm,"Nanq | %.0f | %.0f | ",pospropt[jk],posproptt);
+ fprintf(ficresphtm,"Nanq | %.0f | %.0f | ",pospropt[jk],posproptt);
}else{
- fprintf(ficresphtm,"%.5f | %.0f | %.0f | ",pospropt[jk]/posproptt,pospropt[jk],posproptt);
+ fprintf(ficresphtm,"%.5f | %.0f | %.0f | ",pospropt[jk]/posproptt,pospropt[jk],posproptt);
}
}
fprintf(ficresphtm,"
\n");
@@ -4402,7 +4456,7 @@ Title=%s
Datafile=%s Firstpass=%d La
fprintf(ficresphtmfr,"
\n");
} /* end selected combination of covariate j1 */
dateintmean=dateintsum/k2cpt;
-
+
fclose(ficresp);
fclose(ficresphtm);
fclose(ficresphtmfr);
@@ -4918,7 +4972,7 @@ void concatwav(int wav[], int **dh, int
double ***p3mat;
double eip;
- pstamp(ficreseij);
+ /* pstamp(ficreseij); */
fprintf(ficreseij,"# (a) Life expectancies by health status at initial age and (b) health expectancies by health status at initial age\n");
fprintf(ficreseij,"# Age");
for(i=1; i<=nlstate;i++){
@@ -5283,6 +5337,14 @@ void concatwav(int wav[], int **dh, int
fprintf(ficlog,"Computing total mortality p.j=w1*p1j+w2*p2j+..: result on file '%s' \n",fileresprobmorprev);
pstamp(ficresprobmorprev);
fprintf(ficresprobmorprev,"# probabilities of dying before estepm=%d months for people of exact age and weighted probabilities w1*p1j+w2*p2j+... stand dev in()\n",estepm);
+ fprintf(ficresprobmorprev,"# Selected quantitative variables and dummies");
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ fprintf(ficresprobmorprev," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
+ for(j=1;j<=cptcoveff;j++)
+ fprintf(ficresprobmorprev,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(ij,j)]);
+ fprintf(ficresprobmorprev,"\n");
+
fprintf(ficresprobmorprev,"# Age cov=%-d",ij);
for(j=nlstate+1; j<=(nlstate+ndeath);j++){
fprintf(ficresprobmorprev," p.%-d SE",j);
@@ -5556,7 +5618,9 @@ void concatwav(int wav[], int **dh, int
pstamp(ficresvpl);
fprintf(ficresvpl,"# Standard deviation of period (stable) prevalences \n");
- fprintf(ficresvpl,"# Age");
+ fprintf(ficresvpl,"# Age ");
+ if(nresult >=1)
+ fprintf(ficresvpl," Result# ");
for(i=1; i<=nlstate;i++)
fprintf(ficresvpl," %1d-%1d",i,i);
fprintf(ficresvpl,"\n");
@@ -5642,6 +5706,8 @@ void concatwav(int wav[], int **dh, int
varpl[i][(int)age] = doldm[i][i]; /* Covariances are useless */
fprintf(ficresvpl,"%.0f ",age );
+ if(nresult >=1)
+ fprintf(ficresvpl,"%d ",nres );
for(i=1; i<=nlstate;i++)
fprintf(ficresvpl," %.5f (%.5f)",prlim[i][i],sqrt(varpl[i][(int)age]));
fprintf(ficresvpl,"\n");
@@ -6046,7 +6112,7 @@ void printinghtml(char fileresu[], char
jj1=0;
for(nres=1; nres <= nresult; nres++) /* For each resultline */
- for(k1=1; k1<=m;k1++){
+ for(k1=1; k1<=m;k1++){ /* For each combination of covariate */
if(TKresult[nres]!= k1)
continue;
@@ -6074,51 +6140,51 @@ void printinghtml(char fileresu[], char
}
}
/* aij, bij */
- fprintf(fichtm,"
- Logit model (yours is: 1+age+%s), for example: logit(pij)=log(pij/pii)= aij+ bij age + V1 age + etc. as a function of age: %s_%d-1.svg
\
-",model,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1);
+ fprintf(fichtm,"
- Logit model (yours is: 1+age+%s), for example: logit(pij)=log(pij/pii)= aij+ bij age + V1 age + etc. as a function of age: %s_%d-1-%d.svg
\
+",model,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres);
/* Pij */
- fprintf(fichtm,"
\n- Pij or conditional probabilities to be observed in state j being in state i, %d (stepm) months before: %s_%d-2.svg
\
-",stepm,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1);
+ fprintf(fichtm,"
\n- Pij or conditional probabilities to be observed in state j being in state i, %d (stepm) months before: %s_%d-2-%d.svg
\
+",stepm,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres);
/* Quasi-incidences */
fprintf(fichtm,"
\n- Iij or Conditional probabilities to be observed in state j being in state i %d (stepm) months\
before but expressed in per year i.e. quasi incidences if stepm is small and probabilities too, \
incidence (rates) are the limit when h tends to zero of the ratio of the probability hPij \
-divided by h: hPij/h : %s_%d-3.svg
\
-",stepm,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1);
+divided by h: hPij/h : %s_%d-3-%d.svg
\
+",stepm,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres);
/* Survival functions (period) in state j */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"
\n- Survival functions in state %d. Or probability to survive in state %d being in state (1 to %d) at different ages. %s%d_%d.svg
\
-", cpt, cpt, nlstate, subdirf2(optionfilefiname,"LIJ_"),cpt,jj1,subdirf2(optionfilefiname,"LIJ_"),cpt,jj1,subdirf2(optionfilefiname,"LIJ_"),cpt,jj1);
+ fprintf(fichtm,"
\n- Survival functions in state %d. Or probability to survive in state %d being in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
+", cpt, cpt, nlstate, subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres);
}
/* State specific survival functions (period) */
for(cpt=1; cpt<=nlstate;cpt++){
fprintf(fichtm,"
\n- Survival functions from state %d in each live state and total.\
Or probability to survive in various states (1 to %d) being in state %d at different ages. \
- %s%d_%d.svg
", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,jj1,subdirf2(optionfilefiname,"LIJT_"),cpt,jj1,subdirf2(optionfilefiname,"LIJT_"),cpt,jj1);
+ %s_%d%d-%d.svg
", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres);
}
/* Period (stable) prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"
\n- Convergence to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d.svg
\
-", cpt, cpt, nlstate, subdirf2(optionfilefiname,"P_"),cpt,jj1,subdirf2(optionfilefiname,"P_"),cpt,jj1,subdirf2(optionfilefiname,"P_"),cpt,jj1);
+ fprintf(fichtm,"
\n- Convergence to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
+", cpt, cpt, nlstate, subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres);
}
if(backcast==1){
/* Period (stable) back prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"
\n- Convergence to period (stable) back prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d.svg
\
-", cpt, cpt, nlstate, subdirf2(optionfilefiname,"PB_"),cpt,jj1,subdirf2(optionfilefiname,"PB_"),cpt,jj1,subdirf2(optionfilefiname,"PB_"),cpt,jj1);
+ fprintf(fichtm,"
\n- Convergence to period (stable) back prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
+", cpt, cpt, nlstate, subdirf2(optionfilefiname,"PB_"),cpt,k1,nres,subdirf2(optionfilefiname,"PB_"),cpt,k1,nres,subdirf2(optionfilefiname,"PB_"),cpt,k1,nres);
}
}
if(prevfcast==1){
/* Projection of prevalence up to period (stable) prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"
\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f) up to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s%d_%d.svg
\
-", dateprev1, dateprev2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1);
+ fprintf(fichtm,"
\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f) up to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
+", dateprev1, dateprev2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres);
}
}
for(cpt=1; cpt<=nlstate;cpt++) {
- fprintf(fichtm,"\n
- Life expectancy by health state (%d) at initial age and its decomposition into health expectancies in each alive state (1 to %d) (or area under each survival functions): %s_%d%d.svg
\
-",cpt,nlstate,subdirf2(optionfilefiname,"EXP_"),cpt,jj1,subdirf2(optionfilefiname,"EXP_"),cpt,jj1,subdirf2(optionfilefiname,"EXP_"),cpt,jj1);
+ fprintf(fichtm,"\n
- Life expectancy by health state (%d) at initial age and its decomposition into health expectancies in each alive state (1 to %d) (or area under each survival functions): %s_%d-%d-%d.svg
\
+",cpt,nlstate,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres);
}
/* } /\* end i1 *\/ */
}/* End k1 */
@@ -6177,7 +6243,7 @@ See page 'Matrix of variance-covariance
jj1=0;
- for(nres=1; nres <= nresult; nres++) /* For each resultline */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
for(k1=1; k1<=m;k1++){
if(TKresult[nres]!= k1)
continue;
@@ -6201,17 +6267,18 @@ See page 'Matrix of variance-covariance
}
for(cpt=1; cpt<=nlstate;cpt++) {
fprintf(fichtm,"\n
- Observed (cross-sectional) and period (incidence based) \
-prevalence (with 95%% confidence interval) in state (%d): %s_%d-%d.svg\n
\
-",cpt,subdirf2(optionfilefiname,"V_"),cpt,jj1,subdirf2(optionfilefiname,"V_"),cpt,jj1,subdirf2(optionfilefiname,"V_"),cpt,jj1);
+prevalence (with 95%% confidence interval) in state (%d): %s_%d-%d-%d.svg\n
\
+",cpt,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres);
}
fprintf(fichtm,"\n
- Total life expectancy by age and \
health expectancies in states (1) and (2). If popbased=1 the smooth (due to the model) \
true period expectancies (those weighted with period prevalences are also\
drawn in addition to the population based expectancies computed using\
- observed and cahotic prevalences: %s_%d.svg\n
\
-",subdirf2(optionfilefiname,"E_"),jj1,subdirf2(optionfilefiname,"E_"),jj1,subdirf2(optionfilefiname,"E_"),jj1);
+ observed and cahotic prevalences: %s_%d-%d.svg\n
\
+",subdirf2(optionfilefiname,"E_"),k1,nres,subdirf2(optionfilefiname,"E_"),k1,nres,subdirf2(optionfilefiname,"E_"),k1,nres);
/* } /\* end i1 *\/ */
}/* End k1 */
+ }/* End nres */
fprintf(fichtm,"");
fflush(fichtm);
}
@@ -6270,158 +6337,96 @@ void printinggnuplot(char fileresu[], ch
strcpy(dirfileres,optionfilefiname);
strcpy(optfileres,"vpl");
/* 1eme*/
- for (cpt=1; cpt<= nlstate ; cpt ++) { /* For each live state */
- for (k1=1; k1<= m ; k1 ++) /* For each valid combination of covariate */
+ for (cpt=1; cpt<= nlstate ; cpt ++){ /* For each live state */
+ for (k1=1; k1<= m ; k1 ++){ /* For each valid combination of covariate */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- /* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */
- if(TKresult[nres]!= k1)
- continue;
- /* We are interested in selected combination by the resultline */
- printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);
- fprintf(ficgp,"\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);
- for (k=1; k<=cptcoveff; k++){ /* For each covariate k get corresponding value lv for combination k1 */
- lv= decodtabm(k1,k,cptcoveff); /* Should be the value of the covariate corresponding to k1 combination */
- /* 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]; /* vlv is the value of the covariate lv, 0 or 1 */
- /* For each combination of covariate k1 (V1=1, V3=0), we printed the current covariate k and its value vlv */
- printf(" V%d=%d ",Tvaraff[k],vlv);
- fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
- }
- for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- }
- printf("\n#\n");
- fprintf(ficgp,"\n#\n");
- if(invalidvarcomb[k1]){
- fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
- continue;
- }
+ /* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */
+ if(TKresult[nres]!= k1)
+ continue;
+ /* We are interested in selected combination by the resultline */
+ printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);
+ fprintf(ficgp,"\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);
+ for (k=1; k<=cptcoveff; k++){ /* For each covariate k get corresponding value lv for combination k1 */
+ lv= decodtabm(k1,k,cptcoveff); /* Should be the value of the covariate corresponding to k1 combination */
+ /* 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]; /* vlv is the value of the covariate lv, 0 or 1 */
+ /* For each combination of covariate k1 (V1=1, V3=0), we printed the current covariate k and its value vlv */
+ printf(" V%d=%d ",Tvaraff[k],vlv);
+ fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
+ }
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
+ printf("\n#\n");
+ fprintf(ficgp,"\n#\n");
+ if(invalidvarcomb[k1]){
+ fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
+ continue;
+ }
- fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"V_"),cpt,k1);
- fprintf(ficgp,"\n#set out \"V_%s_%d-%d.svg\" \n",optionfilefiname,cpt,k1);
- fprintf(ficgp,"set xlabel \"Age\" \n\
-set ylabel \"Probability\" \n \
-set ter svg size 640, 480\n \
-plot [%.f:%.f] \"%s\" every :::%d::%d u 1:2 \"%%lf",ageminpar,fage,subdirf2(fileresu,"VPL_"),k1-1,k1-1);
-
- for (i=1; i<= nlstate ; i ++) {
- if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
- else fprintf(ficgp," %%*lf (%%*lf)");
- }
- fprintf(ficgp,"\" t\"Period (stable) prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2+1.96*$3) \"%%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1);
- for (i=1; i<= nlstate ; i ++) {
- if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
- else fprintf(ficgp," %%*lf (%%*lf)");
- }
- fprintf(ficgp,"\" t\"95%% CI\" w l lt 1,\"%s\" every :::%d::%d u 1:($2-1.96*$3) \"%%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1);
- for (i=1; i<= nlstate ; i ++) {
- 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));
- 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 */
- if(cptcoveff ==0){
- fprintf(ficgp,"$%d)) t 'Backward prevalence in state %d' with line ", 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 'Backward prevalence in state %d' ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv], \
- 4+(cpt-1), cpt ); /* 4 or 6 ?*/
- }else{
- fprintf(ficgp,"$%d==%d && $%d==%d && ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv]);
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" \n",subdirf2(optionfilefiname,"V_"),cpt,k1,nres);
+ fprintf(ficgp,"\n#set out \"V_%s_%d-%d-%d.svg\" \n",optionfilefiname,cpt,k1,nres);
+ fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \nset ter svg size 640, 480\nplot [%.f:%.f] \"%s\" every :::%d::%d u 1:($2==%d ? $3:1/0) \"%%lf %%lf",ageminpar,fage,subdirf2(fileresu,"VPL_"),k1-1,k1-1,nres);
+
+ for (i=1; i<= nlstate ; i ++) {
+ if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
+ else fprintf(ficgp," %%*lf (%%*lf)");
+ }
+ fprintf(ficgp,"\" t\"Period (stable) prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2==%d ? $4+1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1,nres);
+ for (i=1; i<= nlstate ; i ++) {
+ if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
+ else fprintf(ficgp," %%*lf (%%*lf)");
+ }
+ fprintf(ficgp,"\" t\"95%% CI\" w l lt 1,\"%s\" every :::%d::%d u 1:($2==%d ? $4-1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1,nres);
+ for (i=1; i<= nlstate ; i ++) {
+ 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));
+ 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 ($2==%d ?$1:1/0):(",subdirf2(fileresu,"PLB_"),nres); /* Age is in 1, nres in 2 to be fixed */
+ if(cptcoveff ==0){
+ fprintf(ficgp,"$%d)) t 'Backward prevalence in state %d' with line ", 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++;
- }
- } /* end covariate */
- } /* end if no covariate */
- } /* end if backcast */
- fprintf(ficgp,"\nset out \n");
+ /* 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 'Backward prevalence in state %d' ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv], \
+ 4+(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 */
+ } /* end if backcast */
+ fprintf(ficgp,"\nset out \n");
+ } /* nres */
} /* k1 */
} /* cpt */
/*2 eme*/
- for (k1=1; k1<= m ; k1 ++)
- for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
- continue;
- fprintf(ficgp,"\n# 2nd: Total life expectancy with CI: 't' files ");
- for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
- lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
- /* 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];
- fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
- }
- for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- }
- fprintf(ficgp,"\n#\n");
- if(invalidvarcomb[k1]){
- fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
- continue;
- }
-
- fprintf(ficgp,"\nset out \"%s_%d.svg\" \n",subdirf2(optionfilefiname,"E_"),k1);
- for(vpopbased=0; vpopbased <= popbased; vpopbased++){ /* Done for vpopbased=0 and vpopbased=1 if popbased==1*/
- if(vpopbased==0)
- fprintf(ficgp,"set ylabel \"Years\" \nset ter svg size 640, 480\nplot [%.f:%.f] ",ageminpar,fage);
- else
- fprintf(ficgp,"\nreplot ");
- for (i=1; i<= nlstate+1 ; i ++) {
- k=2*i;
- fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ?$4 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1, vpopbased);
- for (j=1; j<= nlstate+1 ; j ++) {
- if (j==i) fprintf(ficgp," %%lf (%%lf)");
- else fprintf(ficgp," %%*lf (%%*lf)");
- }
- if (i== 1) fprintf(ficgp,"\" t\"TLE\" w l lt %d, \\\n",i);
- else fprintf(ficgp,"\" t\"LE in state (%d)\" w l lt %d, \\\n",i-1,i+1);
- fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4-$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);
- for (j=1; j<= nlstate+1 ; j ++) {
- if (j==i) fprintf(ficgp," %%lf (%%lf)");
- else fprintf(ficgp," %%*lf (%%*lf)");
- }
- fprintf(ficgp,"\" t\"\" w l lt 0,");
- fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4+$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);
- for (j=1; j<= nlstate+1 ; j ++) {
- if (j==i) fprintf(ficgp," %%lf (%%lf)");
- else fprintf(ficgp," %%*lf (%%*lf)");
- }
- if (i== (nlstate+1)) fprintf(ficgp,"\" t\"\" w l lt 0");
- else fprintf(ficgp,"\" t\"\" w l lt 0,\\\n");
- } /* state */
- } /* vpopbased */
- fprintf(ficgp,"\nset out;set out \"%s_%d.svg\"; replot; set out; \n",subdirf2(optionfilefiname,"E_"),k1); /* Buggy gnuplot */
- } /* k1 end 2 eme*/
-
-
- /*3eme*/
- for (k1=1; k1<= m ; k1 ++)
- for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k)
- continue;
-
- for (cpt=1; cpt<= nlstate ; cpt ++) {
- fprintf(ficgp,"\n# 3d: Life expectancy with EXP_ files: combination=%d state=%d",k1, cpt);
- for (k=1; k<=cptcoveff; k++){ /* For each covariate dummy combination and each value */
+ for (k1=1; k1<= m ; k1 ++){
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
+ fprintf(ficgp,"\n# 2nd: Total life expectancy with CI: 't' files ");
+ for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
/* 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 */
@@ -6431,139 +6436,197 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
}
/* for(k=1; k <= ncovds; k++){ */
for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- }
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
continue;
}
- /* k=2+nlstate*(2*cpt-2); */
- k=2+(nlstate+1)*(cpt-1);
- fprintf(ficgp,"\nset out \"%s_%d%d.svg\" \n",subdirf2(optionfilefiname,"EXP_"),cpt,k1);
- fprintf(ficgp,"set ter svg size 640, 480\n\
+ fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"E_"),k1,nres);
+ for(vpopbased=0; vpopbased <= popbased; vpopbased++){ /* Done for vpopbased=0 and vpopbased=1 if popbased==1*/
+ if(vpopbased==0)
+ fprintf(ficgp,"set ylabel \"Years\" \nset ter svg size 640, 480\nplot [%.f:%.f] ",ageminpar,fage);
+ else
+ fprintf(ficgp,"\nreplot ");
+ for (i=1; i<= nlstate+1 ; i ++) {
+ k=2*i;
+ fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ?$4 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1, vpopbased);
+ for (j=1; j<= nlstate+1 ; j ++) {
+ if (j==i) fprintf(ficgp," %%lf (%%lf)");
+ else fprintf(ficgp," %%*lf (%%*lf)");
+ }
+ if (i== 1) fprintf(ficgp,"\" t\"TLE\" w l lt %d, \\\n",i);
+ else fprintf(ficgp,"\" t\"LE in state (%d)\" w l lt %d, \\\n",i-1,i+1);
+ fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4-$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);
+ for (j=1; j<= nlstate+1 ; j ++) {
+ if (j==i) fprintf(ficgp," %%lf (%%lf)");
+ else fprintf(ficgp," %%*lf (%%*lf)");
+ }
+ fprintf(ficgp,"\" t\"\" w l lt 0,");
+ fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4+$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);
+ for (j=1; j<= nlstate+1 ; j ++) {
+ if (j==i) fprintf(ficgp," %%lf (%%lf)");
+ else fprintf(ficgp," %%*lf (%%*lf)");
+ }
+ if (i== (nlstate+1)) fprintf(ficgp,"\" t\"\" w l lt 0");
+ else fprintf(ficgp,"\" t\"\" w l lt 0,\\\n");
+ } /* state */
+ } /* vpopbased */
+ fprintf(ficgp,"\nset out;set out \"%s_%d.svg\"; replot; set out; \n",subdirf2(optionfilefiname,"E_"),k1); /* Buggy gnuplot */
+ } /* end nres */
+ } /* k1 end 2 eme*/
+
+
+ /*3eme*/
+ for (k1=1; k1<= m ; k1 ++){
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
+
+ for (cpt=1; cpt<= nlstate ; cpt ++) {
+ fprintf(ficgp,"\n# 3d: Life expectancy with EXP_ files: combination=%d state=%d",k1, cpt);
+ for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
+ lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
+ /* 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];
+ fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
+ }
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
+ fprintf(ficgp,"\n#\n");
+ if(invalidvarcomb[k1]){
+ fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
+ continue;
+ }
+
+ /* k=2+nlstate*(2*cpt-2); */
+ k=2+(nlstate+1)*(cpt-1);
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" \n",subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres);
+ fprintf(ficgp,"set ter svg size 640, 480\n\
plot [%.f:%.f] \"%s\" every :::%d::%d u 1:%d t \"e%d1\" w l",ageminpar,fage,subdirf2(fileresu,"E_"),k1-1,k1-1,k,cpt);
- /*fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d-2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);
- for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");
- fprintf(ficgp,"\" t \"e%d1\" w l",cpt);
- fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d+2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);
- for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");
- fprintf(ficgp,"\" t \"e%d1\" w l",cpt);
+ /*fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d-2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);
+ for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");
+ fprintf(ficgp,"\" t \"e%d1\" w l",cpt);
+ fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d+2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);
+ for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");
+ fprintf(ficgp,"\" t \"e%d1\" w l",cpt);
- */
- for (i=1; i< nlstate ; i ++) {
- fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+i,cpt,i+1);
- /* fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileres,"e"),k1-1,k1-1,k+2*i,cpt,i+1);*/
+ */
+ for (i=1; i< nlstate ; i ++) {
+ fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+i,cpt,i+1);
+ /* fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileres,"e"),k1-1,k1-1,k+2*i,cpt,i+1);*/
- }
- fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d.\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+nlstate,cpt);
- }
- }
+ }
+ fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d.\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+nlstate,cpt);
+ }
+ } /* end nres */
+ } /* end kl 3eme */
/* 4eme */
/* Survival functions (period) from state i in state j by initial state i */
- for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
- for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
- continue;
-
- for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
- fprintf(ficgp,"\n#\n#\n# Survival functions in state j : 'LIJ_' files, cov=%d state=%d",k1, cpt);
- for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
- lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
- /* 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];
- fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
- }
- for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- }
- fprintf(ficgp,"\n#\n");
- if(invalidvarcomb[k1]){
- fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
+ for (k1=1; k1<=m; k1++){ /* For each covariate and each value */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
continue;
- }
-
- fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJ_"),cpt,k1);
- fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\
-set ter svg size 640, 480\n \
-unset log y\n \
-plot [%.f:%.f] ", ageminpar, agemaxpar);
- k=3;
- for (i=1; i<= nlstate ; i ++){
- if(i==1){
- fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));
- }else{
- fprintf(ficgp,", '' ");
+ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state cpt*/
+ fprintf(ficgp,"\n#\n#\n# Survival functions in state j : 'LIJ_' files, cov=%d state=%d",k1, cpt);
+ for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
+ lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
+ /* 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];
+ fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
- l=(nlstate+ndeath)*(i-1)+1;
- fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d/($%d",k1,k+l+(cpt-1),k+l);
- for (j=2; j<= nlstate+ndeath ; j ++)
- fprintf(ficgp,"+$%d",k+l+j-1);
- fprintf(ficgp,")) t \"l(%d,%d)\" w l",i,cpt);
- } /* nlstate */
- fprintf(ficgp,"\nset out\n");
- } /* end cpt state*/
- } /* end covariate */
-
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
+ fprintf(ficgp,"\n#\n");
+ if(invalidvarcomb[k1]){
+ fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
+ continue;
+ }
+
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres);
+ fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\
+set ter svg size 640, 480\nunset log y\nplot [%.f:%.f] ", ageminpar, agemaxpar);
+ k=3;
+ for (i=1; i<= nlstate ; i ++){
+ if(i==1){
+ fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));
+ }else{
+ fprintf(ficgp,", '' ");
+ }
+ l=(nlstate+ndeath)*(i-1)+1;
+ fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d/($%d",k1,k+l+(cpt-1),k+l);
+ for (j=2; j<= nlstate+ndeath ; j ++)
+ fprintf(ficgp,"+$%d",k+l+j-1);
+ fprintf(ficgp,")) t \"l(%d,%d)\" w l",i,cpt);
+ } /* nlstate */
+ fprintf(ficgp,"\nset out\n");
+ } /* end cpt state*/
+ } /* end nres */
+ } /* end covariate k1 */
+
/* 5eme */
/* Survival functions (period) from state i in state j by final state j */
- for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
- for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
- continue;
- for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each inital state */
-
- fprintf(ficgp,"\n#\n#\n# Survival functions in state j and all livestates from state i by final state j: 'lij' files, cov=%d state=%d",k1, cpt);
- for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
- lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
- /* 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];
- fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
- }
- for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- }
- fprintf(ficgp,"\n#\n");
- if(invalidvarcomb[k1]){
- fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
+ for (k1=1; k1<= m ; k1++){ /* For each covariate combination if any */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
continue;
- }
+ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each inital state */
+ fprintf(ficgp,"\n#\n#\n# Survival functions in state j and all livestates from state i by final state j: 'lij' files, cov=%d state=%d",k1, cpt);
+ for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
+ lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
+ /* 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];
+ fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
+ }
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
+ fprintf(ficgp,"\n#\n");
+ if(invalidvarcomb[k1]){
+ fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
+ continue;
+ }
- fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJT_"),cpt,k1);
- fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\
-set ter svg size 640, 480\n \
-unset log y\n \
-plot [%.f:%.f] ", ageminpar, agemaxpar);
- k=3;
- for (j=1; j<= nlstate ; j ++){ /* Lived in state j */
- if(j==1)
- fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));
- else
- fprintf(ficgp,", '' ");
- l=(nlstate+ndeath)*(cpt-1) +j;
- fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):($%d",k1,k+l);
- /* for (i=2; i<= nlstate+ndeath ; i ++) */
- /* fprintf(ficgp,"+$%d",k+l+i-1); */
- fprintf(ficgp,") t \"l(%d,%d)\" w l",cpt,j);
- } /* nlstate */
- fprintf(ficgp,", '' ");
- fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):(",k1);
- for (j=1; j<= nlstate ; j ++){ /* Lived in state j */
- l=(nlstate+ndeath)*(cpt-1) +j;
- if(j < nlstate)
- fprintf(ficgp,"$%d +",k+l);
- else
- fprintf(ficgp,"$%d) t\"l(%d,.)\" w l",k+l,cpt);
- }
- fprintf(ficgp,"\nset out\n");
- } /* end cpt state*/
- } /* end covariate */
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres);
+ fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\
+set ter svg size 640, 480\nunset log y\nplot [%.f:%.f] ", ageminpar, agemaxpar);
+ k=3;
+ for (j=1; j<= nlstate ; j ++){ /* Lived in state j */
+ if(j==1)
+ fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));
+ else
+ fprintf(ficgp,", '' ");
+ l=(nlstate+ndeath)*(cpt-1) +j;
+ fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):($%d",k1,k+l);
+ /* for (i=2; i<= nlstate+ndeath ; i ++) */
+ /* fprintf(ficgp,"+$%d",k+l+i-1); */
+ fprintf(ficgp,") t \"l(%d,%d)\" w l",cpt,j);
+ } /* nlstate */
+ fprintf(ficgp,", '' ");
+ fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):(",k1);
+ for (j=1; j<= nlstate ; j ++){ /* Lived in state j */
+ l=(nlstate+ndeath)*(cpt-1) +j;
+ if(j < nlstate)
+ fprintf(ficgp,"$%d +",k+l);
+ else
+ fprintf(ficgp,"$%d) t\"l(%d,.)\" w l",k+l,cpt);
+ }
+ fprintf(ficgp,"\nset out\n");
+ } /* end cpt state*/
+ } /* end covariate */
+ } /* end nres */
/* 6eme */
/* CV preval stable (period) for each covariate */
@@ -6591,11 +6654,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
continue;
}
- fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"P_"),cpt,k1);
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" \n",subdirf2(optionfilefiname,"P_"),cpt,k1,nres);
fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\
-set ter svg size 640, 480\n \
-unset log y\n \
-plot [%.f:%.f] ", ageminpar, agemaxpar);
+set ter svg size 640, 480\nunset log y\nplot [%.f:%.f] ", ageminpar, agemaxpar);
k=3; /* Offset */
for (i=1; i<= nlstate ; i ++){
if(i==1)
@@ -6639,11 +6700,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
continue;
}
- fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"PB_"),cpt,k1);
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" \n",subdirf2(optionfilefiname,"PB_"),cpt,k1,nres);
fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\
-set ter svg size 640, 480\n \
-unset log y\n \
-plot [%.f:%.f] ", ageminpar, agemaxpar);
+set ter svg size 640, 480\nunset log y\nplot [%.f:%.f] ", ageminpar, agemaxpar);
k=3; /* Offset */
for (i=1; i<= nlstate ; i ++){
if(i==1)
@@ -6693,11 +6752,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
}
fprintf(ficgp,"# hpijx=probability over h years, hp.jx is weighted by observed prev\n ");
- fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"PROJ_"),cpt,k1);
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" \n",subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres);
fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Prevalence\" \n\
-set ter svg size 640, 480\n \
-unset log y\n \
-plot [%.f:%.f] ", ageminpar, agemaxpar);
+set ter svg size 640, 480\nunset log y\nplot [%.f:%.f] ", ageminpar, agemaxpar);
for (i=1; i<= nlstate+1 ; i ++){ /* nlstate +1 p11 p21 p.1 */
/*# V1 = 1 V2 = 0 yearproj age p11 p21 p.1 p12 p22 p.2 p13 p23 p.3*/
/*# 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 */
@@ -6764,8 +6821,8 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
} /* End if prevfcast */
- /* proba elementaires */
- fprintf(ficgp,"\n##############\n#MLE estimated parameters\n#############\n");
+ /* 9eme writing MLE parameters */
+ fprintf(ficgp,"\n##############\n#9eme MLE estimated parameters\n#############\n");
for(i=1,jk=1; i <=nlstate; i++){
fprintf(ficgp,"# initial state %d\n",i);
for(k=1; k <=(nlstate+ndeath); k++){
@@ -6782,7 +6839,8 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
fprintf(ficgp,"##############\n#\n");
/*goto avoid;*/
- fprintf(ficgp,"\n##############\n#Graphics of probabilities or incidences\n#############\n");
+ /* 10eme Graphics of probabilities or incidences using written MLE parameters */
+ fprintf(ficgp,"\n##############\n#10eme Graphics of probabilities or incidences\n#############\n");
fprintf(ficgp,"# logi(p12/p11)=a12+b12*age+c12age*age+d12*V1+e12*V1*age\n");
fprintf(ficgp,"# logi(p12/p11)=p1 +p2*age +p3*age*age+ p4*V1+ p5*V1*age\n");
fprintf(ficgp,"# logi(p13/p11)=a13+b13*age+c13age*age+d13*V1+e13*V1*age\n");
@@ -6797,9 +6855,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
fprintf(ficgp,"# +exp(a14+b14*age+c14age*age+d14*V1+e14*V1*age)+...)\n");
fprintf(ficgp,"#\n");
for(ng=1; ng<=3;ng++){ /* Number of graphics: first is logit, 2nd is probabilities, third is incidences per year*/
- fprintf(ficgp,"#Number of graphics: first is logit, 2nd is probabilities, third is incidences per year \n");
+ fprintf(ficgp,"#Number of graphics: first is logit, 2nd is probabilities, third is incidences per year\n");
fprintf(ficgp,"#model=%s \n",model);
- fprintf(ficgp,"# ng=%d\n",ng);
+ fprintf(ficgp,"# Type of graphic ng=%d\n",ng);
fprintf(ficgp,"# jk=1 to 2^%d=%d\n",cptcoveff,m);/* to be checked */
for(jk=1; jk <=m; jk++) /* For each combination of covariate */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
@@ -6810,7 +6868,7 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
}
fprintf(ficgp,"\n#\n");
- fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" ",subdirf2(optionfilefiname,"PE_"),jk,ng);
+ fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" ",subdirf2(optionfilefiname,"PE_"),jk,ng,nres);
fprintf(ficgp,"\nset ter svg size 640, 480 ");
if (ng==1){
fprintf(ficgp,"\nset ylabel \"Value of the logit of the model\"\n"); /* exp(a12+b12*x) could be nice */
@@ -6856,7 +6914,7 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
/* printf("Tage[%d]=%d, j=%d\n", ij, Tage[ij], j); */
if(j==Tage[ij]) { /* Product by age */
if(ij <=cptcovage) { /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, 2 V5 and V1 */
- if(Dummy[j]==0){
+ if(DummyV[j]==0){
fprintf(ficgp,"+p%d*%d*x",i+j+2+nagesqr-1,Tinvresult[nres][Tvar[j]]);;
}else{ /* quantitative */
fprintf(ficgp,"+p%d*%f*x",i+j+2+nagesqr-1,Tqinvresult[nres][Tvar[j]]); /* Tqinvresult in decoderesult */
@@ -6867,8 +6925,8 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
}else if(j==Tprod[ijp]) { /* */
/* printf("Tprod[%d]=%d, j=%d\n", ij, Tprod[ijp], j); */
if(ijp <=cptcovprod) { /* Product */
- if(Dummy[Tvard[ijp][1]]==0){/* Vn is dummy */
- if(Dummy[Tvard[ijp][2]]==0){/* Vn and Vm are dummy */
+ if(DummyV[Tvard[ijp][1]]==0){/* Vn is dummy */
+ if(DummyV[Tvard[ijp][2]]==0){/* Vn and Vm are dummy */
/* fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,nbcode[Tvard[ijp][1]][codtabm(jk,j)],nbcode[Tvard[ijp][2]][codtabm(jk,j)]); */
fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tinvresult[nres][Tvard[ijp][2]]);
}else{ /* Vn is dummy and Vm is quanti */
@@ -6876,12 +6934,13 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);
}
}else{ /* Vn*Vm Vn is quanti */
- if(Dummy[Tvard[ijp][2]]==0){
+ if(DummyV[Tvard[ijp][2]]==0){
fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][2]],Tqinvresult[nres][Tvard[ijp][1]]);
}else{ /* Both quanti */
fprintf(ficgp,"+p%d*%f*%f",i+j+2+nagesqr-1,Tqinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);
}
}
+ ijp++;
}
} else{ /* simple covariate */
/* fprintf(ficgp,"+p%d*%d",i+j+2+nagesqr-1,nbcode[Tvar[j]][codtabm(jk,j)]); /\* Valgrind bug nbcode *\/ */
@@ -7175,8 +7234,7 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
}
for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficresf," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
}
fprintf(ficresf," yearproj age");
for(j=1; j<=nlstate+ndeath;j++){
@@ -7785,14 +7843,36 @@ int readdata(char datafile[], int firsto
/*-------- data file ----------*/
FILE *fic;
char dummy[]=" ";
- int i=0, j=0, n=0, iv=0;
+ int i=0, j=0, n=0, iv=0, v;
int lstra;
int linei, month, year,iout;
char line[MAXLINE], linetmp[MAXLINE];
char stra[MAXLINE], strb[MAXLINE];
char *stratrunc;
+ DummyV=ivector(1,NCOVMAX); /* 1 to 3 */
+ FixedV=ivector(1,NCOVMAX); /* 1 to 3 */
+ for(v=1; v <=ncovcol;v++){
+ DummyV[v]=0;
+ FixedV[v]=0;
+ }
+ for(v=ncovcol+1; v <=ncovcol+nqv;v++){
+ DummyV[v]=1;
+ FixedV[v]=0;
+ }
+ for(v=ncovcol+nqv+1; v <=ncovcol+nqv+ntv;v++){
+ DummyV[v]=0;
+ FixedV[v]=1;
+ }
+ for(v=ncovcol+nqv+ntv+1; v <=ncovcol+nqv+ntv+nqtv;v++){
+ DummyV[v]=1;
+ FixedV[v]=1;
+ }
+ for(v=1; v <=ncovcol+nqv+ntv+nqtv;v++){
+ printf("Covariate type in the data: V%d, DummyV(V%d)=%d, FixedV(V%d)=%d\n",v,v,DummyV[v],v,FixedV[v]);
+ fprintf(ficlog,"Covariate type in the data: V%d, DummyV(V%d)=%d, FixedV(V%d)=%d\n",v,v,DummyV[v],v,FixedV[v]);
+ }
if((fic=fopen(datafile,"r"))==NULL) {
printf("Problem while opening datafile: %s with errno='%s'\n", datafile,strerror(errno));fflush(stdout);
@@ -7822,87 +7902,87 @@ int readdata(char datafile[], int firsto
/* Loops on waves */
for (j=maxwav;j>=1;j--){
for (iv=nqtv;iv>=1;iv--){ /* Loop on time varying quantitative variables */
- cutv(stra, strb, line, ' ');
- if(strb[0]=='.') { /* Missing value */
- lval=-1;
- cotqvar[j][iv][i]=-1; /* 0.0/0.0 */
- cotvar[j][ntv+iv][i]=-1; /* For performance reasons */
- if(isalpha(strb[1])) { /* .m or .d Really Missing value */
- printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value. Exiting.\n", strb, linei,i,line,iv, nqtv, j);
- fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value. Exiting.\n", strb, linei,i,line,iv, nqtv, j);fflush(ficlog);
- return 1;
- }
- }else{
- errno=0;
- /* what_kind_of_number(strb); */
- dval=strtod(strb,&endptr);
- /* if( strb[0]=='\0' || (*endptr != '\0')){ */
- /* if(strb != endptr && *endptr == '\0') */
- /* dval=dlval; */
- /* 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 the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,iv, nqtv, j,maxwav);
- fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line, iv, nqtv, j,maxwav);fflush(ficlog);
- return 1;
- }
- cotqvar[j][iv][i]=dval;
- cotvar[j][ntv+iv][i]=dval;
- }
- strcpy(line,stra);
+ cutv(stra, strb, line, ' ');
+ if(strb[0]=='.') { /* Missing value */
+ lval=-1;
+ cotqvar[j][iv][i]=-1; /* 0.0/0.0 */
+ cotvar[j][ntv+iv][i]=-1; /* For performance reasons */
+ if(isalpha(strb[1])) { /* .m or .d Really Missing value */
+ printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value. Exiting.\n", strb, linei,i,line,iv, nqtv, j);
+ fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value. Exiting.\n", strb, linei,i,line,iv, nqtv, j);fflush(ficlog);
+ return 1;
+ }
+ }else{
+ errno=0;
+ /* what_kind_of_number(strb); */
+ dval=strtod(strb,&endptr);
+ /* if( strb[0]=='\0' || (*endptr != '\0')){ */
+ /* if(strb != endptr && *endptr == '\0') */
+ /* dval=dlval; */
+ /* 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 the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,iv, nqtv, j,maxwav);
+ fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line, iv, nqtv, j,maxwav);fflush(ficlog);
+ return 1;
+ }
+ cotqvar[j][iv][i]=dval;
+ cotvar[j][ntv+iv][i]=dval;
+ }
+ strcpy(line,stra);
}/* end loop ntqv */
for (iv=ntv;iv>=1;iv--){ /* Loop on time varying dummies */
- cutv(stra, strb, line, ' ');
- if(strb[0]=='.') { /* Missing value */
- lval=-1;
- }else{
- 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 the %d th dummy covariate out of %d measured at wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,iv, ntv, j,maxwav);
- fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d dummy covariate out of %d measured wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,iv, ntv,j,maxwav);fflush(ficlog);
- return 1;
- }
- }
- if(lval <-1 || lval >1){
- printf("Error reading data around '%ld' at line number %d for individual %d, '%s'\n \
+ cutv(stra, strb, line, ' ');
+ if(strb[0]=='.') { /* Missing value */
+ lval=-1;
+ }else{
+ 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 the %d th dummy covariate out of %d measured at wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,iv, ntv, j,maxwav);
+ fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d dummy covariate out of %d measured wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,iv, ntv,j,maxwav);fflush(ficlog);
+ return 1;
+ }
+ }
+ if(lval <-1 || lval >1){
+ printf("Error reading data around '%ld' at line number %d for individual %d, '%s'\n \
Should be a value of %d(nth) covariate (0 should be the value for the reference and 1\n \
for the alternative. IMaCh does not build design variables automatically, do it yourself.\n \
- For example, for multinomial values like 1, 2 and 3,\n \
- build V1=0 V2=0 for the reference value (1),\n \
- V1=1 V2=0 for (2) \n \
+ For example, for multinomial values like 1, 2 and 3,\n \
+ build V1=0 V2=0 for the reference value (1),\n \
+ V1=1 V2=0 for (2) \n \
and V1=0 V2=1 for (3). V1=1 V2=1 should not exist and the corresponding\n \
- output of IMaCh is often meaningless.\n \
+ output of IMaCh is often meaningless.\n \
Exiting.\n",lval,linei, i,line,j);
- fprintf(ficlog,"Error reading data around '%ld' at line number %d for individual %d, '%s'\n \
+ fprintf(ficlog,"Error reading data around '%ld' at line number %d for individual %d, '%s'\n \
Should be a value of %d(nth) covariate (0 should be the value for the reference and 1\n \
for the alternative. IMaCh does not build design variables automatically, do it yourself.\n \
- For example, for multinomial values like 1, 2 and 3,\n \
- build V1=0 V2=0 for the reference value (1),\n \
- V1=1 V2=0 for (2) \n \
+ For example, for multinomial values like 1, 2 and 3,\n \
+ build V1=0 V2=0 for the reference value (1),\n \
+ V1=1 V2=0 for (2) \n \
and V1=0 V2=1 for (3). V1=1 V2=1 should not exist and the corresponding\n \
- output of IMaCh is often meaningless.\n \
+ output of IMaCh is often meaningless.\n \
Exiting.\n",lval,linei, i,line,j);fflush(ficlog);
- return 1;
- }
- cotvar[j][iv][i]=(double)(lval);
- strcpy(line,stra);
+ return 1;
+ }
+ cotvar[j][iv][i]=(double)(lval);
+ strcpy(line,stra);
}/* end loop ntv */
/* Statuses at wave */
cutv(stra, strb, line, ' ');
if(strb[0]=='.') { /* Missing value */
- lval=-1;
+ lval=-1;
}else{
- 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);
- return 1;
- }
+ 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);
+ return 1;
+ }
}
s[j][i]=lval;
@@ -8207,7 +8287,7 @@ int decodemodel( char model[], int lasto
* - Tvard[k] p Tvard[1][1]@4 {7, 8, 5, 6} for V7*V8 and V5*V6 .
*/
{
- int i, j, k, ks;
+ int i, j, k, ks, v;
int j1, k1, k2, k3, k4;
char modelsav[80];
char stra[80], strb[80], strc[80], strd[80],stre[80];
@@ -8413,7 +8493,7 @@ 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;k<=cptcovt; k++){ /* or cptocvt */
if (Tvar[k] <=ncovcol && Typevar[k]==0 ){ /* Simple fixed dummy (<=ncovcol) covariates */
Fixed[k]= 0;
@@ -8438,7 +8518,7 @@ Dummy[k] 0=dummy (0 1), 1 quantitative (
TvarFind[ncovf]=k;
TvarFD[ncoveff]=Tvar[k]; /* TvarFD[1]=V1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
TvarFDind[ncoveff]=k; /* TvarFDind[1]=9 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
- }else if( Tvar[k] <=ncovcol+nqv && Typevar[k]==0){ /* Remind that product Vn*Vm are added in k*/ /* Only simple fixed quantitative variable */
+ }else if( Tvar[k] <=ncovcol+nqv && Typevar[k]==0){/* Remind that product Vn*Vm are added in k Only simple fixed quantitative variable */
Fixed[k]= 0;
Dummy[k]= 1;
nqfveff++;
@@ -8491,167 +8571,167 @@ Dummy[k] 0=dummy (0 1), 1 quantitative (
TvarA[ncova]=Tvar[k];
TvarAind[ncova]=k;
if (Tvar[k] <=ncovcol ){ /* Product age with fixed dummy covariatee */
- Fixed[k]= 2;
- Dummy[k]= 2;
- modell[k].maintype= ATYPE;
- modell[k].subtype= APFD;
- /* ncoveff++; */
+ Fixed[k]= 2;
+ Dummy[k]= 2;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APFD;
+ /* ncoveff++; */
}else if( Tvar[k] <=ncovcol+nqv) { /* Remind that product Vn*Vm are added in k*/
- Fixed[k]= 2;
- Dummy[k]= 3;
- modell[k].maintype= ATYPE;
- modell[k].subtype= APFQ; /* Product age * fixed quantitative */
- /* nqfveff++; /\* Only simple fixed quantitative variable *\/ */
+ Fixed[k]= 2;
+ Dummy[k]= 3;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APFQ; /* Product age * fixed quantitative */
+ /* nqfveff++; /\* Only simple fixed quantitative variable *\/ */
}else if( Tvar[k] <=ncovcol+nqv+ntv ){
- Fixed[k]= 3;
- Dummy[k]= 2;
- modell[k].maintype= ATYPE;
- modell[k].subtype= APVD; /* Product age * varying dummy */
- /* ntveff++; /\* Only simple time varying dummy variable *\/ */
+ Fixed[k]= 3;
+ Dummy[k]= 2;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APVD; /* Product age * varying dummy */
+ /* ntveff++; /\* Only simple time varying dummy variable *\/ */
}else if( Tvar[k] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 3;
- Dummy[k]= 3;
- modell[k].maintype= ATYPE;
- modell[k].subtype= APVQ; /* Product age * varying quantitative */
- /* nqtveff++;/\* Only simple time varying quantitative variable *\/ */
+ Fixed[k]= 3;
+ Dummy[k]= 3;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APVQ; /* Product age * varying quantitative */
+ /* nqtveff++;/\* Only simple time varying quantitative variable *\/ */
}
}else if (Typevar[k] == 2) { /* product without age */
k1=Tposprod[k];
if(Tvard[k1][1] <=ncovcol){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 1;
- Dummy[k]= 0;
- modell[k].maintype= FTYPE;
- modell[k].subtype= FPDD; /* Product fixed dummy * fixed dummy */
- ncovf++; /* Fixed variables without age */
- TvarF[ncovf]=Tvar[k];
- TvarFind[ncovf]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv){
- Fixed[k]= 0; /* or 2 ?*/
- Dummy[k]= 1;
- modell[k].maintype= FTYPE;
- modell[k].subtype= FPDQ; /* Product fixed dummy * fixed quantitative */
- ncovf++; /* Varying variables without age */
- TvarF[ncovf]=Tvar[k];
- TvarFind[ncovf]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 0;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDD; /* Product fixed dummy * varying dummy */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDQ; /* Product fixed dummy * varying quantitative */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 1;
+ Dummy[k]= 0;
+ modell[k].maintype= FTYPE;
+ modell[k].subtype= FPDD; /* Product fixed dummy * fixed dummy */
+ ncovf++; /* Fixed variables without age */
+ TvarF[ncovf]=Tvar[k];
+ TvarFind[ncovf]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv){
+ Fixed[k]= 0; /* or 2 ?*/
+ Dummy[k]= 1;
+ modell[k].maintype= FTYPE;
+ modell[k].subtype= FPDQ; /* Product fixed dummy * fixed quantitative */
+ ncovf++; /* Varying variables without age */
+ TvarF[ncovf]=Tvar[k];
+ TvarFind[ncovf]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 0;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDD; /* Product fixed dummy * varying dummy */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product fixed dummy * varying quantitative */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }
}else if(Tvard[k1][1] <=ncovcol+nqv){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 0; /* or 2 ?*/
- Dummy[k]= 1;
- modell[k].maintype= FTYPE;
- modell[k].subtype= FPDQ; /* Product fixed quantitative * fixed dummy */
- ncovf++; /* Fixed variables without age */
- TvarF[ncovf]=Tvar[k];
- TvarFind[ncovf]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDQ; /* Product fixed quantitative * varying dummy */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPQQ; /* Product fixed quantitative * varying quantitative */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 0; /* or 2 ?*/
+ Dummy[k]= 1;
+ modell[k].maintype= FTYPE;
+ modell[k].subtype= FPDQ; /* Product fixed quantitative * fixed dummy */
+ ncovf++; /* Fixed variables without age */
+ TvarF[ncovf]=Tvar[k];
+ TvarFind[ncovf]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product fixed quantitative * varying dummy */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPQQ; /* Product fixed quantitative * varying quantitative */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }
}else if(Tvard[k1][1] <=ncovcol+nqv+ntv){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDD; /* Product time varying dummy * fixed dummy */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDQ; /* Product time varying dummy * fixed quantitative */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 0;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDD; /* Product time varying dummy * time varying dummy */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDQ; /* Product time varying dummy * time varying quantitative */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDD; /* Product time varying dummy * fixed dummy */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying dummy * fixed quantitative */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 0;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDD; /* Product time varying dummy * time varying dummy */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying dummy * time varying quantitative */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }
}else if(Tvard[k1][1] <=ncovcol+nqv+ntv+nqtv){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDQ; /* Product time varying quantitative * fixed dummy */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPQQ; /* Product time varying quantitative * fixed quantitative */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPDQ; /* Product time varying quantitative * time varying dummy */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- modell[k].maintype= VTYPE;
- modell[k].subtype= VPQQ; /* Product time varying quantitative * time varying quantitative */
- ncovv++; /* Varying variables without age */
- TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying quantitative * fixed dummy */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPQQ; /* Product time varying quantitative * fixed quantitative */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying quantitative * time varying dummy */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPQQ; /* Product time varying quantitative * time varying quantitative */
+ ncovv++; /* Varying variables without age */
+ TvarV[ncovv]=Tvar[k];
+ TvarVind[ncovv]=k;
+ }
}else{
- printf("Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
- fprintf(ficlog,"Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
- } /* end k1 */
+ printf("Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
+ fprintf(ficlog,"Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
+ } /*end k1*/
}else{
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]);
@@ -9028,61 +9108,63 @@ int prevalence_limit(double *p, double *
i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */
if (cptcovn < 1){i1=1;}
- for(nres=1; nres <= nresult; nres++) /* For each resultline */
- for(k=1; k<=i1;k++){
- if(TKresult[nres]!= k)
- continue;
+ for(k=1; k<=i1;k++){ /* For each combination k of dummy covariates in the model */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k)
+ continue;
- /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
- /* for(cptcov=1,k=0;cptcov<=1;cptcov++){ */
- //for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){
- /* k=k+1; */
- /* to clean */
- //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
- fprintf(ficrespl,"#******");
- printf("#******");
- fprintf(ficlog,"#******");
- for(j=1;j<=cptcoveff ;j++) {/* all covariates */
- fprintf(ficrespl," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]); /* Here problem for varying dummy*/
- printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- }
- for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
- }
- fprintf(ficrespl,"******\n");
- printf("******\n");
- fprintf(ficlog,"******\n");
- if(invalidvarcomb[k]){
- printf("\nCombination (%d) ignored because no case \n",k);
- fprintf(ficrespl,"#Combination (%d) ignored because no case \n",k);
- fprintf(ficlog,"\nCombination (%d) ignored because no case \n",k);
- continue;
- }
+ /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
+ /* for(cptcov=1,k=0;cptcov<=1;cptcov++){ */
+ //for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){
+ /* k=k+1; */
+ /* to clean */
+ //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
+ fprintf(ficrespl,"#******");
+ printf("#******");
+ fprintf(ficlog,"#******");
+ for(j=1;j<=cptcoveff ;j++) {/* all covariates */
+ fprintf(ficrespl," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]); /* Here problem for varying dummy*/
+ printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ }
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficrespl," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
+ fprintf(ficrespl,"******\n");
+ printf("******\n");
+ fprintf(ficlog,"******\n");
+ if(invalidvarcomb[k]){
+ printf("\nCombination (%d) ignored because no case \n",k);
+ fprintf(ficrespl,"#Combination (%d) ignored because no case \n",k);
+ fprintf(ficlog,"\nCombination (%d) ignored because no case \n",k);
+ continue;
+ }
- fprintf(ficrespl,"#Age ");
- for(j=1;j<=cptcoveff;j++) {
- fprintf(ficrespl,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- }
- for(i=1; i<=nlstate;i++) fprintf(ficrespl," %d-%d ",i,i);
- fprintf(ficrespl,"Total Years_to_converge\n");
+ fprintf(ficrespl,"#Age ");
+ for(j=1;j<=cptcoveff;j++) {
+ fprintf(ficrespl,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ }
+ for(i=1; i<=nlstate;i++) fprintf(ficrespl," %d-%d ",i,i);
+ fprintf(ficrespl,"Total Years_to_converge\n");
- for (age=agebase; age<=agelim; age++){
- /* for (age=agebase; age<=agebase; age++){ */
- prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, ncvyearp, k, nres);
- fprintf(ficrespl,"%.0f ",age );
- for(j=1;j<=cptcoveff;j++)
- fprintf(ficrespl,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- tot=0.;
- for(i=1; i<=nlstate;i++){
- tot += prlim[i][i];
- fprintf(ficrespl," %.5f", prlim[i][i]);
- }
- fprintf(ficrespl," %.3f %d\n", tot, *ncvyearp);
- } /* Age */
- /* was end of cptcod */
- } /* cptcov */
+ for (age=agebase; age<=agelim; age++){
+ /* for (age=agebase; age<=agebase; age++){ */
+ prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, ncvyearp, k, nres);
+ fprintf(ficrespl,"%.0f ",age );
+ for(j=1;j<=cptcoveff;j++)
+ fprintf(ficrespl,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ tot=0.;
+ for(i=1; i<=nlstate;i++){
+ tot += prlim[i][i];
+ fprintf(ficrespl," %.5f", prlim[i][i]);
+ }
+ fprintf(ficrespl," %.3f %d\n", tot, *ncvyearp);
+ } /* Age */
+ /* was end of cptcod */
+ } /* cptcov */
+ } /* nres */
return 0;
}
@@ -9123,69 +9205,69 @@ int back_prevalence_limit(double *p, dou
i1=pow(2,cptcoveff);
if (cptcovn < 1){i1=1;}
- for(nres=1; nres <= nresult; nres++) /* For each resultline */
- for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
- if(TKresult[nres]!= k)
- continue;
- //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
- fprintf(ficresplb,"#******");
- printf("#******");
- fprintf(ficlog,"#******");
- for(j=1;j<=cptcoveff ;j++) {/* all covariates */
- fprintf(ficresplb," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- }
- for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
- printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
- fprintf(ficresplb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
- fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
- }
- fprintf(ficresplb,"******\n");
- printf("******\n");
- fprintf(ficlog,"******\n");
- if(invalidvarcomb[k]){
- printf("\nCombination (%d) ignored because no cases \n",k);
- fprintf(ficresplb,"#Combination (%d) ignored because no cases \n",k);
- fprintf(ficlog,"\nCombination (%d) ignored because no cases \n",k);
- continue;
- }
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
+ if(TKresult[nres]!= k)
+ continue;
+ //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
+ fprintf(ficresplb,"#******");
+ printf("#******");
+ fprintf(ficlog,"#******");
+ for(j=1;j<=cptcoveff ;j++) {/* all covariates */
+ fprintf(ficresplb," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ }
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficresplb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
+ fprintf(ficresplb,"******\n");
+ printf("******\n");
+ fprintf(ficlog,"******\n");
+ if(invalidvarcomb[k]){
+ printf("\nCombination (%d) ignored because no cases \n",k);
+ fprintf(ficresplb,"#Combination (%d) ignored because no cases \n",k);
+ fprintf(ficlog,"\nCombination (%d) ignored because no cases \n",k);
+ continue;
+ }
- fprintf(ficresplb,"#Age ");
- for(j=1;j<=cptcoveff;j++) {
- fprintf(ficresplb,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- }
- for(i=1; i<=nlstate;i++) fprintf(ficresplb," %d-%d ",i,i);
- fprintf(ficresplb,"Total Years_to_converge\n");
+ fprintf(ficresplb,"#Age ");
+ for(j=1;j<=cptcoveff;j++) {
+ fprintf(ficresplb,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ }
+ for(i=1; i<=nlstate;i++) fprintf(ficresplb," %d-%d ",i,i);
+ fprintf(ficresplb,"Total Years_to_converge\n");
- for (age=agebase; age<=agelim; age++){
- /* for (age=agebase; age<=agebase; age++){ */
- if(mobilavproj > 0){
- /* bprevalim(bprlim, mobaverage, nlstate, p, age, ageminpar, agemaxpar, oldm, savm, doldm, dsavm, ftolpl, ncvyearp, k); */
- /* bprevalim(bprlim, mobaverage, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
- bprevalim(bprlim, mobaverage, nlstate, p, age, ftolpl, ncvyearp, k);
- }else if (mobilavproj == 0){
- printf("There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
- fprintf(ficlog,"There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
- exit(1);
- }else{
- /* bprevalim(bprlim, probs, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
- bprevalim(bprlim, probs, nlstate, p, age, ftolpl, ncvyearp, k);
- }
- fprintf(ficresplb,"%.0f ",age );
- for(j=1;j<=cptcoveff;j++)
- fprintf(ficresplb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- tot=0.;
- for(i=1; i<=nlstate;i++){
- tot += bprlim[i][i];
- fprintf(ficresplb," %.5f", bprlim[i][i]);
- }
- fprintf(ficresplb," %.3f %d\n", tot, *ncvyearp);
- } /* Age */
- /* was end of cptcod */
- } /* cptcov */
-
+ for (age=agebase; age<=agelim; age++){
+ /* for (age=agebase; age<=agebase; age++){ */
+ if(mobilavproj > 0){
+ /* bprevalim(bprlim, mobaverage, nlstate, p, age, ageminpar, agemaxpar, oldm, savm, doldm, dsavm, ftolpl, ncvyearp, k); */
+ /* bprevalim(bprlim, mobaverage, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
+ bprevalim(bprlim, mobaverage, nlstate, p, age, ftolpl, ncvyearp, k);
+ }else if (mobilavproj == 0){
+ printf("There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
+ fprintf(ficlog,"There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
+ exit(1);
+ }else{
+ /* bprevalim(bprlim, probs, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
+ bprevalim(bprlim, probs, nlstate, p, age, ftolpl, ncvyearp, k);
+ }
+ fprintf(ficresplb,"%.0f ",age );
+ for(j=1;j<=cptcoveff;j++)
+ fprintf(ficresplb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ tot=0.;
+ for(i=1; i<=nlstate;i++){
+ tot += bprlim[i][i];
+ fprintf(ficresplb," %.5f", bprlim[i][i]);
+ }
+ fprintf(ficresplb," %.3f %d\n", tot, *ncvyearp);
+ } /* Age */
+ /* was end of cptcod */
+ } /* end of any combination */
+ } /* end of nres */
/* hBijx(p, bage, fage); */
/* fclose(ficrespijb); */
@@ -9277,7 +9359,7 @@ int hPijx(double *p, int bage, int fage)
int ageminl;
int hstepm;
int nhstepm;
- int h, i, i1, j, k;
+ int h, i, i1, j, k, nres;
double agedeb;
double ***p3mat;
@@ -9305,48 +9387,54 @@ int hPijx(double *p, int bage, int fage)
/* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
/* /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */
/* k=k+1; */
- for (k=1; k <= (int) pow(2,cptcoveff); k++){
- fprintf(ficrespijb,"\n#****** ");
- for(j=1;j<=cptcoveff;j++)
- fprintf(ficrespijb,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- fprintf(ficrespijb,"******\n");
- if(invalidvarcomb[k]){
- fprintf(ficrespijb,"\n#Combination (%d) ignored because no cases \n",k);
- continue;
- }
-
- /* for (agedeb=fage; agedeb>=bage; agedeb--){ /\* If stepm=6 months *\/ */
- for (agedeb=bage; agedeb<=fage; agedeb++){ /* If stepm=6 months and estepm=24 (2 years) */
- /* nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm); /\* Typically 20 years = 20*12/6=40 *\/ */
- nhstepm=(int) rint((agedeb-ageminl)*YEARM/stepm); /* Typically 20 years = 20*12/6=40 */
- nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 */
-
- /* nhstepm=nhstepm*YEARM; aff par mois*/
-
- p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
- /* oldm=oldms;savm=savms; */
- /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k); */
- hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm, k);
- /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm,oldm,savm, dnewm, doldm, dsavm, k); */
- fprintf(ficrespijb,"# Cov Agex agex-h hpijx with i,j=");
- for(i=1; i<=nlstate;i++)
- for(j=1; j<=nlstate+ndeath;j++)
- fprintf(ficrespijb," %1d-%1d",i,j);
- fprintf(ficrespijb,"\n");
- for (h=0; h<=nhstepm; h++){
- /*agedebphstep = agedeb + h*hstepm/YEARM*stepm;*/
- fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb - h*hstepm/YEARM*stepm );
- /* fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb + h*hstepm/YEARM*stepm ); */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
+ if(TKresult[nres]!= k)
+ continue;
+ fprintf(ficrespijb,"\n#****** ");
+ for(j=1;j<=cptcoveff;j++)
+ fprintf(ficrespijb,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ fprintf(ficrespijb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
+ fprintf(ficrespijb,"******\n");
+ if(invalidvarcomb[k]){
+ fprintf(ficrespijb,"\n#Combination (%d) ignored because no cases \n",k);
+ continue;
+ }
+
+ /* for (agedeb=fage; agedeb>=bage; agedeb--){ /\* If stepm=6 months *\/ */
+ for (agedeb=bage; agedeb<=fage; agedeb++){ /* If stepm=6 months and estepm=24 (2 years) */
+ /* nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm); /\* Typically 20 years = 20*12/6=40 *\/ */
+ nhstepm=(int) rint((agedeb-ageminl)*YEARM/stepm); /* Typically 20 years = 20*12/6=40 */
+ nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 */
+
+ /* nhstepm=nhstepm*YEARM; aff par mois*/
+
+ p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ /* oldm=oldms;savm=savms; */
+ /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k); */
+ hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm, k);
+ /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm,oldm,savm, dnewm, doldm, dsavm, k); */
+ fprintf(ficrespijb,"# Cov Agex agex-h hpijx with i,j=");
for(i=1; i<=nlstate;i++)
for(j=1; j<=nlstate+ndeath;j++)
- fprintf(ficrespijb," %.5f", p3mat[i][j][h]);
+ fprintf(ficrespijb," %1d-%1d",i,j);
fprintf(ficrespijb,"\n");
- }
- free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
- fprintf(ficrespijb,"\n");
- }
- /*}*/
- }
+ for (h=0; h<=nhstepm; h++){
+ /*agedebphstep = agedeb + h*hstepm/YEARM*stepm;*/
+ fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb - h*hstepm/YEARM*stepm );
+ /* fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb + h*hstepm/YEARM*stepm ); */
+ for(i=1; i<=nlstate;i++)
+ for(j=1; j<=nlstate+ndeath;j++)
+ fprintf(ficrespijb," %.5f", p3mat[i][j][h]);
+ fprintf(ficrespijb,"\n");
+ }
+ free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ fprintf(ficrespijb,"\n");
+ } /* end age deb */
+ } /* end combination */
+ } /* end nres */
return 0;
} /* hBijx */
@@ -10751,6 +10839,7 @@ Please run with mle=-1 to get a correct
ungetc(c,ficpar);
fgets(line, MAXLINE, ficpar);
fputs(line,stdout);
+ fputs(line,ficres);
fputs(line,ficparo);
}
ungetc(c,ficpar);
@@ -10767,6 +10856,7 @@ Please run with mle=-1 to get a correct
fgets(line, MAXLINE, ficpar);
fputs(line,stdout);
fputs(line,ficparo);
+ fputs(line,ficres);
}
ungetc(c,ficpar);
@@ -10785,14 +10875,17 @@ Please run with mle=-1 to get a correct
fputs(line,stdout);
fputs(line,ficparo);
fputs(line,ficlog);
+ fputs(line,ficres);
continue;
}else
break;
}
+ if (!feof(ficpar))
while((num_filled=sscanf(line,"result:%[^\n]\n",resultline)) !=EOF){
- if (num_filled == 0)
+ if (num_filled == 0){
resultline[0]='\0';
- else if (num_filled != 1){
+ break;
+ } else if (num_filled != 1){
printf("ERROR %d: result line should be at minimum 'result=' %s\n",num_filled, line);
}
nresult++; /* Sum of resultlines */
@@ -10803,12 +10896,16 @@ Please run with mle=-1 to get a correct
goto end;
}
decoderesult(resultline, nresult); /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */
+ fprintf(ficparo,"result: %s\n",resultline);
+ fprintf(ficres,"result: %s\n",resultline);
+ fprintf(ficlog,"result: %s\n",resultline);
while(fgets(line, MAXLINE, ficpar)) {
/* If line starts with a # it is a comment */
if (line[0] == '#') {
numlinepar++;
fputs(line,stdout);
fputs(line,ficparo);
+ fputs(line,ficres);
fputs(line,ficlog);
continue;
}else
@@ -10818,7 +10915,7 @@ Please run with mle=-1 to get a correct
break;
else{ /* Processess output results for this combination of covariate values */
}
- }
+ } /* end while */
@@ -10953,6 +11050,8 @@ Please run with mle=-1 to get a correct
}
printf("Computing Health Expectancies: result on file '%s' ...", filerese);fflush(stdout);
fprintf(ficlog,"Computing Health Expectancies: result on file '%s' ...", filerese);fflush(ficlog);
+
+ pstamp(ficreseij);
i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */
if (cptcovn < 1){i1=1;}
@@ -11065,6 +11164,7 @@ Please run with mle=-1 to get a correct
fprintf(ficrescveij,"******\n");
fprintf(ficresvij,"\n#****** ");
+ /* pstamp(ficresvij); */
for(j=1;j<=cptcoveff;j++)
fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
@@ -11242,6 +11342,8 @@ Please run with mle=-1 to get a correct
free_ivector(ncodemaxwundef,1,NCOVMAX);
free_ivector(Dummy,-1,NCOVMAX);
free_ivector(Fixed,-1,NCOVMAX);
+ free_ivector(DummyV,1,NCOVMAX);
+ free_ivector(FixedV,1,NCOVMAX);
free_ivector(Typevar,-1,NCOVMAX);
free_ivector(Tvar,1,NCOVMAX);
free_ivector(TvarsQ,1,NCOVMAX);