--- imach/src/imach.c 2022/08/31 08:23:16 1.335 +++ imach/src/imach.c 2022/09/04 17:40:33 1.338 @@ -1,6 +1,23 @@ -/* $Id: imach.c,v 1.335 2022/08/31 08:23:16 brouard Exp $ +/* $Id: imach.c,v 1.338 2022/09/04 17:40:33 brouard Exp $ $State: Exp $ $Log: imach.c,v $ + Revision 1.338 2022/09/04 17:40:33 brouard + Summary: 0.99r36 + + * imach.c (Module): Now the easy runs i.e. without result or + model=1+age only did not work. The defautl combination should be 1 + and not 0 because everything hasn't been tranformed yet. + + Revision 1.337 2022/09/02 14:26:02 brouard + Summary: version 0.99r35 + + * src/imach.c: Version 0.99r35 because it outputs same results with + 1+age+V1+V1*age for females and 1+age for females only + (education=1 noweight) + + Revision 1.336 2022/08/31 09:52:36 brouard + *** empty log message *** + Revision 1.335 2022/08/31 08:23:16 brouard Summary: improvements... @@ -1284,12 +1301,12 @@ typedef struct { #define ODIRSEPARATOR '\\' #endif -/* $Id: imach.c,v 1.335 2022/08/31 08:23:16 brouard Exp $ */ +/* $Id: imach.c,v 1.338 2022/09/04 17:40:33 brouard Exp $ */ /* $State: Exp $ */ #include "version.h" char version[]=__IMACH_VERSION__; -char copyright[]="August 2022,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2020, Nihon University 2021-202, INED 2000-2022"; -char fullversion[]="$Revision: 1.335 $ $Date: 2022/08/31 08:23:16 $"; +char copyright[]="September 2022,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2020, Nihon University 2021-202, INED 2000-2022"; +char fullversion[]="$Revision: 1.338 $ $Date: 2022/09/04 17:40:33 $"; char strstart[80]; char optionfilext[10], optionfilefiname[FILENAMELENGTH]; int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */ @@ -1505,8 +1522,9 @@ int **nbcode, *Tvar; /**< model=V2 => Tv /*Tvar[k]= 5 4 3 6 5 2 7 1 1 */ /* nsd 1 2 3 */ /* Counting single dummies covar fixed or tv */ /*TnsdVar[Tvar] 1 2 3 */ +/*Tvaraff[nsd] 4 3 1 */ /* ID of single dummy cova fixed or timevary*/ /*TvarsD[nsd] 4 3 1 */ /* ID of single dummy cova fixed or timevary*/ -/*TvarsDind[k] 2 3 9 */ /* position K of single dummy cova */ +/*TvarsDind[nsd] 2 3 9 */ /* position K of single dummy cova */ /* nsq 1 2 */ /* Counting single quantit tv */ /* TvarsQ[k] 5 2 */ /* Number of single quantitative cova */ /* TvarsQind 1 6 */ /* position K of single quantitative cova */ @@ -2499,7 +2517,8 @@ void powell(double p[], double **xi, int xits=vector(1,n); *fret=(*func)(p); for (j=1;j<=n;j++) pt[j]=p[j]; - rcurr_time = time(NULL); + rcurr_time = time(NULL); + fp=(*fret); /* Initialisation */ for (*iter=1;;++(*iter)) { ibig=0; del=0.0; @@ -2507,8 +2526,10 @@ void powell(double p[], double **xi, int /* (void) gettimeofday(&curr_time,&tzp); */ rcurr_time = time(NULL); curr_time = *localtime(&rcurr_time); - printf("\nPowell iter=%d -2*LL=%.12f gain=%.12f=%.3g %ld sec. %ld sec.",*iter,*fret, fp-*fret,fp-*fret, rcurr_time-rlast_time, rcurr_time-rstart_time);fflush(stdout); - fprintf(ficlog,"\nPowell iter=%d -2*LL=%.12f gain=%.12f=%.3g %ld sec. %ld sec.",*iter,*fret, fp-*fret,fp-*fret,rcurr_time-rlast_time, rcurr_time-rstart_time); fflush(ficlog); + /* printf("\nPowell iter=%d -2*LL=%.12f gain=%.12f=%.3g %ld sec. %ld sec.",*iter,*fret, fp-*fret,fp-*fret, rcurr_time-rlast_time, rcurr_time-rstart_time);fflush(stdout); */ + /* fprintf(ficlog,"\nPowell iter=%d -2*LL=%.12f gain=%.12f=%.3g %ld sec. %ld sec.",*iter,*fret, fp-*fret,fp-*fret,rcurr_time-rlast_time, rcurr_time-rstart_time); fflush(ficlog); */ + printf("\nPowell iter=%d -2*LL=%.12f gain=%.3lg %ld sec. %ld sec.",*iter,*fret,fp-*fret, rcurr_time-rlast_time, rcurr_time-rstart_time);fflush(stdout); + fprintf(ficlog,"\nPowell iter=%d -2*LL=%.12f gain=%.3lg %ld sec. %ld sec.",*iter,*fret,fp-*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); */ fp=(*fret); /* From former iteration or initial value */ for (i=1;i<=n;i++) { @@ -2811,7 +2832,7 @@ void powell(double p[], double **xi, int double **prevalim(double **prlim, int nlstate, double x[], double age, double **oldm, double **savm, double ftolpl, int *ncvyear, int ij, int nres) { - /**< Computes the prevalence limit in each live state at age x and for covariate combination ij + /**< Computes the prevalence limit in each live state at age x and for covariate combination ij . Nicely done * (and selected quantitative values in nres) * by left multiplying the unit * matrix by transitions matrix until convergence is reached with precision ftolpl @@ -3508,7 +3529,8 @@ double **matprod2(double **out, double * double ***hpxij(double ***po, int nhstepm, double age, int hstepm, double *x, int nlstate, int stepm, double **oldm, double **savm, int ij, int nres ) { - /* Computes the transition matrix starting at age 'age' and dummies values in each resultline (loop on ij to find the corresponding combination) to over + /* Already optimized with precov. + Computes the transition matrix starting at age 'age' and dummies values in each resultline (loop on ij to find the corresponding combination) to over 'nhstepm*hstepm*stepm' months (i.e. until age (in years) age+nhstepm*hstepm*stepm/12) by multiplying nhstepm*hstepm matrices. @@ -3839,16 +3861,17 @@ double ***hbxij(double ***po, int nhstep /*************** log-likelihood *************/ double func( double *x) { - int i, ii, j, k, mi, d, kk; + int i, ii, j, k, mi, d, kk, kf=0; int ioffset=0; double l, ll[NLSTATEMAX+1], cov[NCOVMAX+1]; double **out; double lli; /* Individual log likelihood */ int s1, s2; int iv=0, iqv=0, itv=0, iqtv=0 ; /* Index of varying covariate, fixed quantitative cov, time varying covariate, quantitative time varying covariate */ + double bbh, survp; - long ipmx; double agexact; + double agebegin, ageend; /*extern weight */ /* We are differentiating ll according to initial status */ /* for (i=1;i<=npar;i++) printf("%f ", x[i]);*/ @@ -3871,12 +3894,12 @@ double func( double *x) */ ioffset=2+nagesqr ; /* Fixed */ - for (k=1; k<=ncovf;k++){ /* For each fixed covariate dummu or quant or prod */ + for (kf=1; kf<=ncovf;kf++){ /* For each fixed covariate dummu or quant or prod */ /* # V1=sex, V2=raedyrs Quant Fixed, State=livarnb4..livarnb11, V3=iadl4..iald11, V4=adlw4..adlw11, V5=r4bmi..r11bmi */ /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* TvarF[1]=Tvar[6]=2, TvarF[2]=Tvar[7]=7, TvarF[3]=Tvar[9]=1 ID of fixed covariates or product V2, V1*V2, V1 */ /* TvarFind; TvarFind[1]=6, TvarFind[2]=7, TvarFind[3]=9 *//* Inverse V2(6) is first fixed (single or prod) */ - cov[ioffset+TvarFind[k]]=covar[Tvar[TvarFind[k]]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V1 is fixed (TvarFind[1]=6)*/ + cov[ioffset+TvarFind[kf]]=covar[Tvar[TvarFind[kf]]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V1 is fixed (TvarFind[1]=6)*/ /* V1*V2 (7) TvarFind[2]=7, TvarFind[3]=9 */ } /* In model V2+V1*V4+age*V3+V3*V2 Tvar[1] is V2, Tvar[2=V1*V4] @@ -3891,7 +3914,8 @@ double func( double *x) But if the variable is not in the model TTvar[iv] is the real variable effective in the model: meaning that decodemodel should be used cotvar[mw[mi+1][i]][TTvar[iv]][i] */ - for(mi=1; mi<= wav[i]-1; mi++){ + for(mi=1; mi<= wav[i]-1; mi++){ /* Varying with waves */ + /* Wave varying (but not age varying) */ for(k=1; k <= ncovv ; k++){ /* Varying covariates in the model (single and product but no age )"V5+V4+V3+V4*V3+V5*age+V1*age+V1" +TvarVind 1,2,3,4(V4*V3) Tvar[1]@7{5, 4, 3, 6, 5, 1, 1 ; 6 because the created covar is after V5 and is 6, minus 1+1, 3,2,1,4 positions in cotvar*/ /* cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]][i]; but where is the crossproduct? */ cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]-ncovcol-nqv][i]; @@ -3901,6 +3925,9 @@ double func( double *x) oldm[ii][j]=(ii==j ? 1.0 : 0.0); savm[ii][j]=(ii==j ? 1.0 : 0.0); } + + agebegin=agev[mw[mi][i]][i]; /* Age at beginning of effective wave */ + ageend=agev[mw[mi][i]][i] + (dh[mi][i])*stepm/YEARM; /* Age at end of effective wave and at the end of transition */ for(d=0; d(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*(savm[s1][s2])):log((1.+bbh)*out[s1][s2]));*/ /* linear interpolation */ @@ -4198,6 +4225,11 @@ double funcone( double *x) for(k=1; k<=nlstate; k++) ll[k]=0.; ioffset=0; for (i=1,ipmx=0, sw=0.; i<=imx; i++){ + /* Computes the values of the ncovmodel covariates of the model + depending if the covariates are fixed or varying (age dependent) and stores them in cov[] + Then computes with function pmij which return a matrix p[i][j] giving the elementary probability + to be observed in j being in i according to the model. + */ /* ioffset=2+nagesqr+cptcovage; */ ioffset=2+nagesqr; /* Fixed */ @@ -4215,6 +4247,19 @@ double funcone( double *x) /* cov[2+9]=covar[Tvar[9]][i]; */ /* cov[2+9]=covar[1][i]; V1 */ } + /* In model V2+V1*V4+age*V3+V3*V2 Tvar[1] is V2, Tvar[2=V1*V4] + is 5, Tvar[3=age*V3] should not be computed because of age Tvar[4=V3*V2]=6 + has been calculated etc */ + /* For an individual i, wav[i] gives the number of effective waves */ + /* We compute the contribution to Likelihood of each effective transition + mw[mi][i] is real wave of the mi th effectve wave */ + /* Then statuses are computed at each begin and end of an effective wave s1=s[ mw[mi][i] ][i]; + s2=s[mw[mi+1][i]][i]; + And the iv th varying covariate is the cotvar[mw[mi+1][i]][iv][i] + But if the variable is not in the model TTvar[iv] is the real variable effective in the model: + meaning that decodemodel should be used cotvar[mw[mi+1][i]][TTvar[iv]][i] + */ + /* This part may be useless now because everythin should be in covar */ /* for (k=1; k<=nqfveff;k++){ /\* Simple and product fixed Quantitative covariates without age* products *\/ */ /* cov[++ioffset]=coqvar[TvarFQ[k]][i];/\* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V2 and V1*V2 is fixed (k=6 and 7?)*\/ */ /* } */ @@ -4272,7 +4317,19 @@ double funcone( double *x) savm=oldm; oldm=newm; } /* end mult */ - + /*lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/ /* Original formula */ + /* But now since version 0.9 we anticipate for bias at large stepm. + * If stepm is larger than one month (smallest stepm) and if the exact delay + * (in months) between two waves is not a multiple of stepm, we rounded to + * the nearest (and in case of equal distance, to the lowest) interval but now + * we keep into memory the bias bh[mi][i] and also the previous matrix product + * (i.e to dh[mi][i]-1) saved in 'savm'. Then we inter(extra)polate the + * probability in order to take into account the bias as a fraction of the way + * from savm to out if bh is negative or even beyond if bh is positive. bh varies + * -stepm/2 to stepm/2 . + * For stepm=1 the results are the same as for previous versions of Imach. + * For stepm > 1 the results are less biased than in previous versions. + */ s1=s[mw[mi][i]][i]; s2=s[mw[mi+1][i]][i]; /* if(s2==-1){ */ @@ -5080,7 +5137,7 @@ Title=%s
Datafile=%s Firstpass=%d La /* }else */ /* TODO TODO codtabm(j1,z1) or codtabm(j1,Tvaraff[z1]]z1)*/ /* if( iind >=imx-3) printf("Searching error iind=%d Tvaraff[z1]=%d covar[Tvaraff[z1]][iind]=%.f TnsdVar[Tvaraff[z1]]=%d, cptcoveff=%d, cptcovs=%d \n",iind, Tvaraff[z1], covar[Tvaraff[z1]][iind],TnsdVar[Tvaraff[z1]],cptcoveff, cptcovs); */ if(Tvaraff[z1]<1 || Tvaraff[z1]>=NCOVMAX) - printf("Error Tvaraff[z1]=%d<1 or >=%d, cptcoveff=%d model=%s\n",Tvaraff[z1],NCOVMAX, cptcoveff, model); + printf("Error Tvaraff[z1]=%d<1 or >=%d, cptcoveff=%d model=1+age+%s\n",Tvaraff[z1],NCOVMAX, cptcoveff, model); if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]){ /* for combination j1 of covariates */ /* Tests if the value of the covariate z1 for this individual iind responded to combination j1 (V4=1 V3=0) */ bool=0; /* bool should be equal to 1 to be selected, one covariate value failed */ @@ -5390,7 +5447,7 @@ Title=%s
Datafile=%s Firstpass=%d La printf("# This combination (%d) is not valid and no result will be produced\n",j1); invalidvarcomb[j1]=1; }else{ - fprintf(ficresphtm,"\n

This combination (%d) is valid and result will be produced.

",j1); + fprintf(ficresphtm,"\n

This combination (%d) is valid and result will be produced (or no resultline).

",j1); invalidvarcomb[j1]=0; } fprintf(ficresphtmfr,"\n"); @@ -6234,6 +6291,7 @@ void concatwav(int wav[], int **dh, int /* Covariances of health expectancies eij and of total life expectancies according to initial status i, ei. . */ + /* Very time consuming function, but already optimized with precov */ int i, j, nhstepm, hstepm, h, nstepm, k, cptj, cptj2, i2, j2, ij, ji; int nhstepma, nstepma; /* Decreasing with age */ double age, agelim, hf; @@ -6502,11 +6560,16 @@ void concatwav(int wav[], int **dh, int 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 */ /* To be done*/ - fprintf(ficresprobmorprev," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); + + /* We use TinvDoQresult[nres][resultmodel[nres][j] we sort according to the equation model and the resultline: it is a choice */ + /* for (j=1; j<= nsq; j++){ /\* For each selected (single) quantitative value *\/ /\* To be done*\/ */ + /* fprintf(ficresprobmorprev," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* } */ + for (j=1; j<= cptcovs; j++){ /* For each selected (single) quantitative value */ /* To be done*/ + fprintf(ficresprobmorprev," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); } - for(j=1;j<=cptcoveff;j++) - fprintf(ficresprobmorprev," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(ij,TnsdVar[Tvaraff[j]])]); + /* for(j=1;j<=cptcoveff;j++) */ + /* fprintf(ficresprobmorprev," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(ij,TnsdVar[Tvaraff[j]])]); */ fprintf(ficresprobmorprev,"\n"); fprintf(ficresprobmorprev,"# Age cov=%-d",ij); @@ -7153,9 +7216,9 @@ To be simple, these graphs help to under /* Including quantitative variables of the resultline to be done */ for (z1=1; z1<=cptcovs; z1++){ /* Loop on each variable of this resultline */ - printf("Varprob modelresult[%d][%d]=%d model=%s \n",nres, z1, modelresult[nres][z1], model); - fprintf(ficlog,"Varprob modelresult[%d][%d]=%d model=%s \n",nres, z1, modelresult[nres][z1], model); - /* fprintf(ficlog,"Varprob modelresult[%d][%d]=%d model=%s resultline[%d]=%s \n",nres, z1, modelresult[nres][z1], model, nres, resultline[nres]); */ + printf("Varprob modelresult[%d][%d]=%d model=1+age+%s \n",nres, z1, modelresult[nres][z1], model); + fprintf(ficlog,"Varprob modelresult[%d][%d]=%d model=1+age+%s \n",nres, z1, modelresult[nres][z1], model); + /* fprintf(ficlog,"Varprob modelresult[%d][%d]=%d model=1+age+%s resultline[%d]=%s \n",nres, z1, modelresult[nres][z1], model, nres, resultline[nres]); */ if(Dummy[modelresult[nres][z1]]==0){/* Dummy variable of the variable in position modelresult in the model corresponding to z1 in resultline */ if(Fixed[modelresult[nres][z1]]==0){ /* Fixed referenced to model equation */ fprintf(ficresprob,"V%d=%d ",Tvresult[nres][z1],Tresult[nres][z1]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ @@ -7177,7 +7240,7 @@ To be simple, these graphs help to under } }else if(Dummy[modelresult[nres][z1]]==1){ /* Quanti variable */ /* For each selected (single) quantitative value */ - fprintf(ficresprob," V%d=%f ",Tvqresult[nres][z1],Tqresult[nres][z1]); + fprintf(ficresprob," V%d=%lg ",Tvqresult[nres][z1],Tqresult[nres][z1]); if(Fixed[modelresult[nres][z1]]==0){ /* Fixed */ fprintf(ficresprob,"fixed "); fprintf(ficresprobcov,"fixed "); @@ -7524,29 +7587,39 @@ void printinghtml(char fileresu[], char jj1=0; fprintf(fichtm," \n"); fprintf(fichtm,"\ @@ -7708,7 +7780,7 @@ See page 'Matrix of variance-covariance /* - Population forecasting (if popforecast=1): pop%s
\n */ /*
",fileres,fileres,fileres,fileres); */ /* else */ -/* fprintf(fichtm,"\n No population forecast: popforecast = %d (instead of 1) or stepm = %d (instead of 1) or model=%s (instead of .)

\n",popforecast, stepm, model); */ +/* fprintf(fichtm,"\n No population forecast: popforecast = %d (instead of 1) or stepm = %d (instead of 1) or model=1+age+%s (instead of .)

\n",popforecast, stepm, model); */ fflush(fichtm); m=pow(2,cptcoveff); @@ -7719,28 +7791,24 @@ See page 'Matrix of variance-covariance jj1=0; fprintf(fichtm," \n"); jj1=0; for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - for(k1=1; k1<=m;k1++){ - if(m != 1 && TKresult[nres]!= k1) - continue; + /* k1=nres; */ + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* for(k1=1; k1<=m;k1++){ */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ /* for(i1=1; i1<=ncodemax[k1];i1++){ */ jj1++; if (cptcovn > 0) { fprintf(fichtm,"\n

"); fprintf(fichtm,"

************ Results for covariates"); - for (cpt=1; cpt<=cptcoveff;cpt++){ /**< cptcoveff number of variables */ - fprintf(fichtm," V%d=%d ",Tvresult[nres][cpt],Tresult[nres][cpt]); - printf(" V%d=%d ",Tvresult[nres][cpt],Tresult[nres][cpt]);fflush(stdout); + for (cpt=1; cpt<=cptcovs;cpt++){ /**< cptcoveff number of variables */ + fprintf(fichtm," V%d=%lg ",Tvresult[nres][cpt],TinvDoQresult[nres][Tvresult[nres][cpt]]); + printf(" V%d=%lg ",Tvresult[nres][cpt],TinvDoQresult[nres][Tvresult[nres][cpt]]); /* fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]); */ } - for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */ - fprintf(fichtm," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); - } - fprintf(fichtm," (model=%s) ************\n
",model); + fprintf(fichtm," (model=1+age+%s) ************\n
",model); if(invalidvarcomb[k1]){ fprintf(fichtm,"\n

Combination (%d) ignored because no cases

\n",k1); continue; } - } + } /* If cptcovn >0 */ for(cpt=1; cpt<=nlstate;cpt++) { fprintf(fichtm,"\n
- Observed (cross-sectional with mov_average=%d) and period (incidence based) \ prevalence (with 95%% confidence interval) in state (%d): %s_%d-%d-%d.svg",mobilav,cpt,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres); @@ -7800,7 +7865,6 @@ true period expectancies (those weighted fprintf(fichtm," (data from text file %s.txt) \n
",subdirf2(optionfilefiname,"T_"),subdirf2(optionfilefiname,"T_")); fprintf(fichtm,"",subdirf2(optionfilefiname,"E_"),k1,nres); /* } /\* end i1 *\/ */ - }/* End k1 */ }/* End nres */ fprintf(fichtm,""); fflush(fichtm); @@ -7878,31 +7942,37 @@ void printinggnuplot(char fileresu[], ch 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 (k1=1; k1<= m ; k1 ++){ /\* For each valid combination of covariate *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ /* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */ - if(m != 1 && TKresult[nres]!= k1) - continue; + /* if(m != 1 && 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: Forward (stable period) prevalence with CI: 'VPL_' files and live state =%d ", cpt); strcpy(gplotlabel,"("); - 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 *\/ */ - lv=codtabm(k1,TnsdVar[Tvaraff[k]]); - /* 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); - sprintf(gplotlabel+strlen(gplotlabel)," 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]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + + /* 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 *\\/ *\/ */ + /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */ + /* /\* 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); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); /* printf("\n#\n"); */ @@ -7916,7 +7986,7 @@ void printinggnuplot(char fileresu[], ch 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 label \"Alive state %d %s\" at graph 0.98,0.5 center rotate font \"Helvetica,12\"\n",cpt,gplotlabel); */ - fprintf(ficgp,"set title \"Alive state %d %s model=%s\" font \"Helvetica,12\"\n",cpt,gplotlabel,model); + fprintf(ficgp,"set title \"Alive state %d %s model=1+age+%s\" font \"Helvetica,12\"\n",cpt,gplotlabel,model); 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_"),nres-1,nres-1,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); */ /* k1-1 error should be nres-1*/ @@ -8015,33 +8085,38 @@ void printinggnuplot(char fileresu[], ch /* fprintf(ficgp,"\nset out ;unset label;\n"); */ fprintf(ficgp,"\nset out ;unset title;\n"); } /* nres */ - } /* k1 */ + /* } /\* k1 *\/ */ } /* cpt */ /*2 eme*/ - for (k1=1; k1<= m ; k1 ++){ + /* for (k1=1; k1<= m ; k1 ++){ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ fprintf(ficgp,"\n# 2nd: Total life expectancy with CI: 't' files "); strcpy(gplotlabel,"("); - 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 *\/ */ - lv=codtabm(k1,TnsdVar[Tvaraff[k]]); - /* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%d ",Tvaraff[k],vlv); - } - /* 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]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* 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 *\\/ *\/ */ + /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */ + /* /\* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%d ",Tvaraff[k],vlv); */ + /* } */ + /* /\* 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]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); @@ -8083,33 +8158,38 @@ void printinggnuplot(char fileresu[], ch } /* vpopbased */ fprintf(ficgp,"\nset out;set out \"%s_%d-%d.svg\"; replot; set out; unset label;\n",subdirf2(optionfilefiname,"E_"),k1,nres); /* Buggy gnuplot */ } /* end nres */ - } /* k1 end 2 eme*/ + /* } /\* k1 end 2 eme*\/ */ /*3eme*/ - for (k1=1; k1<= m ; k1 ++){ + /* for (k1=1; k1<= m ; k1 ++){ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ for (cpt=1; cpt<= nlstate ; cpt ++) { /* Fragile no verification of covariate values */ fprintf(ficgp,"\n\n# 3d: Life expectancy with EXP_ files: combination=%d state=%d",k1, cpt); strcpy(gplotlabel,"("); - 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 *\/ */ - lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /* Should be the covariate value corresponding to combination k1 and covariate k */ - /* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][resultmodel[nres][k4]]); - } + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* 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 *\\/ *\/ */ + /* lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /\* Should be the covariate value corresponding to combination k1 and covariate k *\/ */ + /* /\* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][resultmodel[nres][k4]]); */ + } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); if(invalidvarcomb[k1]){ @@ -8140,31 +8220,36 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u } fprintf(ficgp,"\nunset label;\n"); } /* end nres */ - } /* end kl 3eme */ + /* } /\* 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 and each value */ + /* for (k1=1; k1<=m; k1++){ /\* For each covariate and each value *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state cpt*/ strcpy(gplotlabel,"("); - 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=codtabm(k1,TnsdVar[Tvaraff[k]]); - /* 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]; */ - vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); + fprintf(ficgp,"\n#\n#\n# Survival functions in state %d : 'LIJ_' files, cov=%d state=%d", cpt, k1, cpt); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* for (k=1; k<=cptcoveff; k++){ /\* For each covariate and each value *\/ */ + /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */ + /* /\* 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]; *\/ */ + /* vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); @@ -8193,31 +8278,36 @@ set ter svg size 640, 480\nunset log y\n fprintf(ficgp,"\nset out; unset label;\n"); } /* end cpt state*/ } /* end nres */ - } /* end covariate k1 */ + /* } /\* 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 (k1=1; k1<= m ; k1++){ /\* For each covariate combination if any *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each inital state */ strcpy(gplotlabel,"("); 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=codtabm(k1,TnsdVar[Tvaraff[k]]); - /* 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]; */ - vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* for (k=1; k<=cptcoveff; k++){ /\* For each covariate and each value *\/ */ + /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */ + /* /\* 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]; *\/ */ + /* vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); @@ -8253,32 +8343,37 @@ set ter svg size 640, 480\nunset log y\n } fprintf(ficgp,"\nset out; unset label;\n"); } /* end cpt state*/ - } /* end covariate */ + /* } /\* end covariate *\/ */ } /* end nres */ /* 6eme */ /* CV preval stable (period) for each covariate */ - for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ + /* for (k1=1; k1<= m ; k1 ++) /\* For each covariate combination if any *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state of arrival */ strcpy(gplotlabel,"("); fprintf(ficgp,"\n#\n#\n#CV preval stable (forward): 'pij' files, covariatecombination#=%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 *\/ */ - lv=codtabm(k1,TnsdVar[Tvaraff[k]]); - /* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* 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 *\\/ *\/ */ + /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */ + /* /\* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); @@ -8311,27 +8406,32 @@ set ter svg size 640, 480\nunset log y\n /* 7eme */ if(prevbcast == 1){ /* CV backward prevalence for each covariate */ - for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ + /* for (k1=1; k1<= m ; k1 ++) /\* For each covariate combination if any *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life origin state */ strcpy(gplotlabel,"("); fprintf(ficgp,"\n#\n#\n#CV Backward stable prevalence: 'pijb' files, covariatecombination#=%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 *\/ */ - lv=codtabm(k1,TnsdVar[Tvaraff[k]]); - /* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* 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 *\\/ *\/ */ + /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */ + /* /\* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); @@ -8369,27 +8469,32 @@ set ter svg size 640, 480\nunset log y\n if(prevfcast==1){ /* Projection from cross-sectional to forward stable (period) prevalence for each covariate */ - for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ + /* for (k1=1; k1<= m ; k1 ++) /\* For each covariate combination if any *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */ strcpy(gplotlabel,"("); fprintf(ficgp,"\n#\n#\n#Projection of prevalence to forward stable prevalence (period): 'PROJ_' files, covariatecombination#=%d state=%d",k1, cpt); - for (k=1; k<=cptcoveff; k++){ /* For each correspondig covariate value */ - /* lv= decodtabm(k1,k,cptcoveff); /\* Should be the covariate value corresponding to k1 combination and kth covariate *\/ */ - lv=codtabm(k1,TnsdVar[Tvaraff[k]]); - /* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* for (k=1; k<=cptcoveff; k++){ /\* For each correspondig covariate value *\/ */ + /* /\* lv= decodtabm(k1,k,cptcoveff); /\\* Should be the covariate value corresponding to k1 combination and kth covariate *\\/ *\/ */ + /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */ + /* /\* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); @@ -8486,27 +8591,32 @@ set ter svg size 640, 480\nunset log y\n if(prevbcast==1){ /* Back projection from cross-sectional to stable (mixed) for each covariate */ - for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ + /* for (k1=1; k1<= m ; k1 ++) /\* For each covariate combination if any *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */ strcpy(gplotlabel,"("); fprintf(ficgp,"\n#\n#\n#Back projection of prevalence to stable (mixed) back prevalence: 'BPROJ_' files, covariatecombination#=%d originstate=%d",k1, cpt); - for (k=1; k<=cptcoveff; k++){ /* For each correspondig covariate value */ - /* lv= decodtabm(k1,k,cptcoveff); /\* Should be the covariate value corresponding to k1 combination and kth covariate *\/ */ - lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /* Should be the covariate value corresponding to combination k1 and covariate k */ - /* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); + for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + /* for (k=1; k<=cptcoveff; k++){ /\* For each correspondig covariate value *\/ */ + /* /\* lv= decodtabm(k1,k,cptcoveff); /\\* Should be the covariate value corresponding to k1 combination and kth covariate *\\/ *\/ */ + /* lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /\* Should be the covariate value corresponding to combination k1 and covariate k *\/ */ + /* /\* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ } strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); @@ -8561,19 +8671,25 @@ set ter svg size 640, 480\nunset log y\n fprintf(ficgp," u %d:(",ioffset); kl=0; strcpy(gplotcondition,"("); - for (k=1; k<=cptcoveff; k++){ /* For each covariate writing the chain of conditions */ - /* lv= decodtabm(k1,k,cptcoveff); /\* Should be the covariate value corresponding to combination k1 and covariate k *\/ */ - lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /* Should be the covariate value corresponding to combination k1 and covariate k */ - /* 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]; /\* Value of the modality of Tvaraff[k] *\/ */ - vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - kl++; - sprintf(gplotcondition+strlen(gplotcondition),"$%d==%d && $%d==%d " ,kl,Tvaraff[k], kl+1, nbcode[Tvaraff[k]][lv]); - kl++; - if(k 1) - sprintf(gplotcondition+strlen(gplotcondition)," && "); + for (k=1; k<=cptcovs; k++){ /* For each covariate k of the resultline, get corresponding value lv for combination k1 */ + if(Dummy[modelresult[nres][k]]==0){ /* To be verified */ + /* for (k=1; k<=cptcoveff; k++){ /\* For each covariate writing the chain of conditions *\/ */ + /* lv= decodtabm(k1,k,cptcoveff); /\* Should be the covariate value corresponding to combination k1 and covariate k *\/ */ + /* lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /\* Should be the covariate value corresponding to combination k1 and covariate k *\/ */ + lv=Tvresult[nres][k]; + vlv=TinvDoQresult[nres][Tvresult[nres][k]]; + /* 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]; /\* Value of the modality of Tvaraff[k] *\/ */ + /* vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + kl++; + /* sprintf(gplotcondition+strlen(gplotcondition),"$%d==%d && $%d==%d " ,kl,Tvaraff[k], kl+1, nbcode[Tvaraff[k]][lv]); */ + sprintf(gplotcondition+strlen(gplotcondition),"$%d==%d && $%d==%lg " ,kl,Tvresult[nres][k], kl+1,TinvDoQresult[nres][Tvresult[nres][k]]); + kl++; + if(k 1) + sprintf(gplotcondition+strlen(gplotcondition)," && "); + } } strcpy(gplotcondition+strlen(gplotcondition),")"); /* 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 *\/ */ @@ -8638,31 +8754,43 @@ set ter svg size 640, 480\nunset log y\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,"#model=%s \n",model); + fprintf(ficgp,"#model=1+age+%s \n",model); fprintf(ficgp,"# Type of graphic ng=%d\n",ng); fprintf(ficgp,"# k1=1 to 2^%d=%d\n",cptcoveff,m);/* to be checked */ - for(k1=1; k1 <=m; k1++) /* For each combination of covariate */ + /* for(k1=1; k1 <=m; k1++) /\* For each combination of covariate *\/ */ for(nres=1; nres <= nresult; nres++){ /* For each resultline */ - if(m != 1 && TKresult[nres]!= k1) - continue; - fprintf(ficgp,"\n\n# Combination of dummy k1=%d which is ",k1); + /* k1=nres; */ + k1=TKresult[nres]; + if(TKresult[nres]==0) k1=1; /* To be checked for noresult */ + fprintf(ficgp,"\n\n# Resultline k1=%d ",k1); strcpy(gplotlabel,"("); /*sprintf(gplotlabel+strlen(gplotlabel)," Dummy combination %d ",k1);*/ - for (k=1; k<=cptcoveff; k++){ /* For each correspondig covariate value */ - /* lv= decodtabm(k1,k,cptcoveff); /\* Should be the covariate value corresponding to k1 combination and kth covariate *\/ */ - lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /* Should be the covariate value corresponding to combination k1 and covariate k */ - /* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; - fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); - sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); - } + for (k=1; k<=cptcovs; k++){ /**< cptcovs number of SIMPLE covariates in the model V2+V1 =2 (dummy or quantit or time varying) */ + /* for each resultline nres, and position k, Tvresult[nres][k] gives the name of the variable and + TinvDoQresult[nres][Tvresult[nres][k]] gives its value double or integer) */ + fprintf(ficgp," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + sprintf(gplotlabel+strlen(gplotlabel)," V%d=%lg ",Tvresult[nres][k],TinvDoQresult[nres][Tvresult[nres][k]]); + } + /* if(m != 1 && TKresult[nres]!= k1) */ + /* continue; */ + /* fprintf(ficgp,"\n\n# Combination of dummy k1=%d which is ",k1); */ + /* strcpy(gplotlabel,"("); */ + /* /\*sprintf(gplotlabel+strlen(gplotlabel)," Dummy combination %d ",k1);*\/ */ + /* for (k=1; k<=cptcoveff; k++){ /\* For each correspondig covariate value *\/ */ + /* /\* lv= decodtabm(k1,k,cptcoveff); /\\* Should be the covariate value corresponding to k1 combination and kth covariate *\\/ *\/ */ + /* lv= codtabm(k1,TnsdVar[Tvaraff[k]]); /\* Should be the covariate value corresponding to combination k1 and covariate k *\/ */ + /* /\* 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= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */ + /* fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," 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][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* sprintf(gplotlabel+strlen(gplotlabel)," V%d=%f ",Tvqresult[nres][resultmodel[nres][k4]],Tqresult[nres][resultmodel[nres][k4]]); */ + /* } */ strcpy(gplotlabel+strlen(gplotlabel),")"); fprintf(ficgp,"\n#\n"); fprintf(ficgp,"\nset out \"%s_%d-%d-%d.svg\" ",subdirf2(optionfilefiname,"PE_"),k1,ng,nres); @@ -8878,7 +9006,7 @@ set ter svg size 640, 480\nunset log y\n } /* end k2 */ /* fprintf(ficgp,"\n set out; unset label;set key default;\n"); */ fprintf(ficgp,"\n set out; unset title;set key default;\n"); - } /* end k1 */ + } /* end resultline */ } /* end ng */ /* avoid: */ fflush(ficgp); @@ -9447,23 +9575,27 @@ void prevforecast(char fileres[], double i1=pow(2,cptcoveff); if (cptcovn < 1){i1=1;} - for(nres=1; nres <= nresult; nres++) /* For each resultline */ - for(k=1; k<=i1;k++){ /* We find the combination equivalent to result line values of dummies */ + for(nres=1; nres <= nresult; nres++){ /* For each resultline */ + k=TKresult[nres]; + if(TKresult[nres]==0) k=1; /* To be checked for noresult */ + /* for(k=1; k<=i1;k++){ /\* We find the combination equivalent to result line values of dummies *\/ */ if(i1 != 1 && TKresult[nres]!= k) continue; fprintf(ficresvpl,"\n#****** "); printf("\n#****** "); fprintf(ficlog,"\n#****** "); - for(j=1;j<=cptcoveff;j++) { - fprintf(ficresvpl,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - fprintf(ficlog,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); + for(j=1;j<=cptcovs;j++) { + fprintf(ficresvpl,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + fprintf(ficlog,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + printf("V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + /* fprintf(ficlog,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ } - for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */ - printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficresvpl," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - } + /* for (j=1; j<= nsq; j++){ /\* For each selected (single) quantitative value *\/ */ + /* printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* fprintf(ficresvpl," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* } */ fprintf(ficresvpl,"******\n"); printf("******\n"); fprintf(ficlog,"******\n"); @@ -9504,22 +9636,28 @@ void prevforecast(char fileres[], double i1=pow(2,cptcoveff); if (cptcovn < 1){i1=1;} - for(nres=1; nres <= nresult; nres++) /* For each resultline */ - for(k=1; k<=i1;k++){ - if(i1 != 1 && TKresult[nres]!= k) - continue; + for(nres=1; nres <= nresult; nres++){ /* For each resultline */ + k=TKresult[nres]; + if(TKresult[nres]==0) k=1; /* To be checked for noresult */ + /* for(k=1; k<=i1;k++){ */ + /* if(i1 != 1 && TKresult[nres]!= k) */ + /* continue; */ fprintf(ficresvbl,"\n#****** "); printf("\n#****** "); fprintf(ficlog,"\n#****** "); - for(j=1;j<=cptcoveff;j++) { - fprintf(ficresvbl,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - fprintf(ficlog,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - } - for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */ - printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficresvbl," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); + for (j=1; j<= cptcovs; j++){ /* For each selected (single) quantitative value */ + printf(" V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][resultmodel[nres][j]]); + fprintf(ficresvbl," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][resultmodel[nres][j]]); + fprintf(ficlog," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][resultmodel[nres][j]]); + /* for(j=1;j<=cptcoveff;j++) { */ + /* fprintf(ficresvbl,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* fprintf(ficlog,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* } */ + /* for (j=1; j<= nsq; j++){ /\* For each selected (single) quantitative value *\/ */ + /* printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* fprintf(ficresvbl," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ } fprintf(ficresvbl,"******\n"); printf("******\n"); @@ -10406,8 +10544,8 @@ int decoderesult( char resultline[], int } } if(match == 0){ - printf("Error in result line (Dummy single): V%d is missing in result: %s according to model=%s. Tvar[k1=%d]=%d is different from Tvarsel[k2=%d]=%d.\n",Tvar[k1], resultline, model,k1, Tvar[k1], k2, Tvarsel[k2]); - fprintf(ficlog,"Error in result line (Dummy single): V%d is missing in result: %s according to model=%s\n",Tvar[k1], resultline, model); + printf("Error in result line (Dummy single): V%d is missing in result: %s according to model=1+age+%s. Tvar[k1=%d]=%d is different from Tvarsel[k2=%d]=%d.\n",Tvar[k1], resultline, model,k1, Tvar[k1], k2, Tvarsel[k2]); + fprintf(ficlog,"Error in result line (Dummy single): V%d is missing in result: %s according to model=1+age+%s\n",Tvar[k1], resultline, model); return 1; } }else if(Typevar[k1]==1){ /* Product with age We want to get the position k2 in the resultline of the product k1 in the model line*/ @@ -10423,8 +10561,8 @@ int decoderesult( char resultline[], int } } if(match == 0){ - printf("Error in result line (Product with age): V%d is missing in result: %s according to model=%s (Tvarsel[k2=%d]=%d)\n",Tvar[k1], resultline, model, k2, Tvarsel[k2]); - fprintf(ficlog,"Error in result line (Product with age): V%d is missing in result: %s according to model=%s (Tvarsel[k2=%d]=%d)\n",Tvar[k1], resultline, model, k2, Tvarsel[k2]); + printf("Error in result line (Product with age): V%d is missing in result: %s according to model=1+age+%s (Tvarsel[k2=%d]=%d)\n",Tvar[k1], resultline, model, k2, Tvarsel[k2]); + fprintf(ficlog,"Error in result line (Product with age): V%d is missing in result: %s according to model=1+age+%s (Tvarsel[k2=%d]=%d)\n",Tvar[k1], resultline, model, k2, Tvarsel[k2]); return 1; } }else if(Typevar[k1]==2){ /* Product No age We want to get the position in the resultline of the product in the model line*/ @@ -10439,8 +10577,8 @@ int decoderesult( char resultline[], int } } if(match == 0){ - printf("Error in result line (Product without age first variable): V%d is missing in result: %s according to model=%s\n",Tvardk[k1][1], resultline, model); - fprintf(ficlog,"Error in result line (Product without age first variable): V%d is missing in result: %s according to model=%s\n",Tvardk[k1][1], resultline, model); + printf("Error in result line (Product without age first variable): V%d is missing in result: %s according to model=1+age+%s\n",Tvardk[k1][1], resultline, model); + fprintf(ficlog,"Error in result line (Product without age first variable): V%d is missing in result: %s according to model=1+age+%s\n",Tvardk[k1][1], resultline, model); return 1; } match=0; @@ -10453,8 +10591,8 @@ int decoderesult( char resultline[], int } } if(match == 0){ - printf("Error in result line (Product without age second variable): V%d is missing in result: %s according to model=%s\n",Tvardk[k1][2], resultline, model); - fprintf(ficlog,"Error in result line (Product without age second variable): V%d is missing in result : %s according to model=%s\n",Tvardk[k1][2], resultline, model); + printf("Error in result line (Product without age second variable): V%d is missing in result: %s according to model=1+age+%s\n",Tvardk[k1][2], resultline, model); + fprintf(ficlog,"Error in result line (Product without age second variable): V%d is missing in result : %s according to model=1+age+%s\n",Tvardk[k1][2], resultline, model); return 1; } }/* End of testing */ @@ -10474,12 +10612,12 @@ int decoderesult( char resultline[], int } } if(match == 0){ - printf("Error in result line: variable V%d is missing in model; result: %s, model=%s\n",Tvarsel[k2], resultline, model); - fprintf(ficlog,"Error in result line: variable V%d is missing in model; result: %s, model=%s\n",Tvarsel[k2], resultline, model); + printf("Error in result line: variable V%d is missing in model; result: %s, model=1+age+%s\n",Tvarsel[k2], resultline, model); + fprintf(ficlog,"Error in result line: variable V%d is missing in model; result: %s, model=1+age+%s\n",Tvarsel[k2], resultline, model); return 1; }else if(match > 1){ - printf("Error in result line: %d doubled; result: %s, model=%s\n",k2, resultline, model); - fprintf(ficlog,"Error in result line: %d doubled; result: %s, model=%s\n",k2, resultline, model); + printf("Error in result line: %d doubled; result: %s, model=1+age+%s\n",k2, resultline, model); + fprintf(ficlog,"Error in result line: %d doubled; result: %s, model=1+age+%s\n",k2, resultline, model); return 1; } } @@ -10611,18 +10749,18 @@ int decodemodel( char model[], int lasto return 1; } if (strstr(model,"v") !=0){ - printf("Error. 'v' must be in upper case 'V' model=%s ",model); - fprintf(ficlog,"Error. 'v' must be in upper case model=%s ",model);fflush(ficlog); + printf("Error. 'v' must be in upper case 'V' model=1+age+%s ",model); + fprintf(ficlog,"Error. 'v' must be in upper case model=1+age+%s ",model);fflush(ficlog); return 1; } strcpy(modelsav,model); if ((strpt=strstr(model,"age*age")) !=0){ - printf(" strpt=%s, model=%s\n",strpt, model); + printf(" strpt=%s, model=1+age+%s\n",strpt, model); if(strpt != model){ - printf("Error in model: 'model=%s'; 'age*age' should in first place before other covariates\n \ + printf("Error in model: 'model=1+age+%s'; 'age*age' should in first place before other covariates\n \ 'model=1+age+age*age+V1.' or 'model=1+age+age*age+V1+V1*age.', please swap as well as \n \ corresponding column of parameters.\n",model); - fprintf(ficlog,"Error in model: 'model=%s'; 'age*age' should in first place before other covariates\n \ + fprintf(ficlog,"Error in model: 'model=1+age+%s'; 'age*age' should in first place before other covariates\n \ 'model=1+age+age*age+V1.' or 'model=1+age+age*age+V1+V1*age.', please swap as well as \n \ corresponding column of parameters.\n",model); fflush(ficlog); return 1; @@ -11071,16 +11209,16 @@ Dummy[k] 0=dummy (0 1), 1 quantitative ( if((Typevar[k1]==Typevar[k2]) && (Fixed[k1]==Fixed[k2]) && (Dummy[k1]==Dummy[k2] )){ if((Typevar[k1] == 0 || Typevar[k1] == 1)){ /* Simple or age product */ if(Tvar[k1]==Tvar[k2]){ - printf("Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[k1],Dummy[k1]); - fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[k1],Dummy[k1]); fflush(ficlog); + printf("Error duplication in the model=1+age+%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[k1],Dummy[k1]); + fprintf(ficlog,"Error duplication in the model=1+age+%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[k1],Dummy[k1]); fflush(ficlog); return(1); } }else if (Typevar[k1] ==2){ k3=Tposprod[k1]; k4=Tposprod[k2]; if( ((Tvard[k3][1]== Tvard[k4][1])&&(Tvard[k3][2]== Tvard[k4][2])) || ((Tvard[k3][1]== Tvard[k4][2])&&(Tvard[k3][2]== Tvard[k4][1])) ){ - printf("Error duplication in the model=%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); - fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); fflush(ficlog); + printf("Error duplication in the model=1+age+%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); + fprintf(ficlog,"Error duplication in the model=1+age+%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); fflush(ficlog); return(1); } } @@ -11436,10 +11574,12 @@ int prevalence_limit(double *p, double * i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */ if (cptcovn < 1){i1=1;} - for(k=1; k<=i1;k++){ /* For each combination k of dummy covariates in the model */ + /* 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(i1 != 1 && TKresult[nres]!= k) - continue; + k=TKresult[nres]; + if(TKresult[nres]==0) k=1; /* To be checked for noresult */ + /* if(i1 != 1 && TKresult[nres]!= k) /\* We found the combination k corresponding to the resultline value of dummies *\/ */ + /* continue; */ /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */ /* for(cptcov=1,k=0;cptcov<=1;cptcov++){ */ @@ -11450,17 +11590,19 @@ int prevalence_limit(double *p, double * fprintf(ficrespl,"#******"); printf("#******"); fprintf(ficlog,"#******"); - for(j=1;j<=cptcoveff ;j++) {/* all covariates */ + for(j=1;j<=cptcovs ;j++) {/**< cptcovs number of SIMPLE covariates in the model or resultline V2+V1 =2 (dummy or quantit or time varying) */ /* fprintf(ficrespl," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,Tvaraff[j])]); /\* Here problem for varying dummy*\/ */ - fprintf(ficrespl," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); /* Here problem for varying dummy*/ - printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[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]); - } + /* printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + fprintf(ficrespl," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + printf(" V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + fprintf(ficlog," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][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"); @@ -11472,18 +11614,24 @@ int prevalence_limit(double *p, double * } fprintf(ficrespl,"#Age "); - for(j=1;j<=cptcoveff;j++) { - fprintf(ficrespl,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); + /* for(j=1;j<=cptcoveff;j++) { */ + /* fprintf(ficrespl,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* } */ + for(j=1;j<=cptcovs;j++) { /* New the quanti variable is added */ + fprintf(ficrespl,"V%d %lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][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); + /**< Computes the prevalence limit in each live state at age x and for covariate combination (k and) nres */ + prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, ncvyearp, k, nres); /* Nicely done */ fprintf(ficrespl,"%.0f ",age ); - for(j=1;j<=cptcoveff;j++) - fprintf(ficrespl,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); + /* for(j=1;j<=cptcoveff;j++) */ + /* fprintf(ficrespl,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + for(j=1;j<=cptcovs;j++) + fprintf(ficrespl,"%d %lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); tot=0.; for(i=1; i<=nlstate;i++){ tot += prlim[i][i]; @@ -11492,8 +11640,8 @@ int prevalence_limit(double *p, double * fprintf(ficrespl," %.3f %d\n", tot, *ncvyearp); } /* Age */ /* was end of cptcod */ - } /* cptcov */ - } /* nres */ + } /* nres */ + /* } /\* for each combination *\/ */ return 0; } @@ -11535,23 +11683,30 @@ int back_prevalence_limit(double *p, dou 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(i1 != 1 && TKresult[nres]!= k) - continue; - /*printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));*/ + /* for(k=1; k<=i1;k++){ /\* For any combination of dummy covariates, fixed and varying *\/ */ + k=TKresult[nres]; + if(TKresult[nres]==0) k=1; /* To be checked for noresult */ + /* if(i1 != 1 && 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,TnsdVar[Tvaraff[j]])]); - printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - } - for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */ - printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficresplb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - } + for(j=1;j<=cptcovs ;j++) {/**< cptcovs number of SIMPLE covariates in the model or resultline V2+V1 =2 (dummy or quantit or time varying) */ + printf(" V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + fprintf(ficresplb," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + fprintf(ficlog," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + } + /* for(j=1;j<=cptcoveff ;j++) {/\* all covariates *\/ */ + /* fprintf(ficresplb," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* } */ + /* for (j=1; j<= nsq; j++){ /\* For each selected (single) quantitative value *\/ */ + /* printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* fprintf(ficresplb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* } */ fprintf(ficresplb,"******\n"); printf("******\n"); fprintf(ficlog,"******\n"); @@ -11563,8 +11718,8 @@ int back_prevalence_limit(double *p, dou } fprintf(ficresplb,"#Age "); - for(j=1;j<=cptcoveff;j++) { - fprintf(ficresplb,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); + for(j=1;j<=cptcovs;j++) { + fprintf(ficresplb,"V%d %lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); } for(i=1; i<=nlstate;i++) fprintf(ficresplb," %d-%d ",i,i); fprintf(ficresplb,"Total Years_to_converge\n"); @@ -11587,8 +11742,8 @@ int back_prevalence_limit(double *p, dou /* exit(1); */ } fprintf(ficresplb,"%.0f ",age ); - for(j=1;j<=cptcoveff;j++) - fprintf(ficresplb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); + for(j=1;j<=cptcovs;j++) + fprintf(ficresplb,"%d %lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); tot=0.; for(i=1; i<=nlstate;i++){ tot += bprlim[i][i]; @@ -11598,7 +11753,7 @@ int back_prevalence_limit(double *p, dou } /* Age */ /* was end of cptcod */ /*fprintf(ficresplb,"\n");*/ /* Seems to be necessary for gnuplot only if two result lines and no covariate. */ - } /* end of any combination */ + /* } /\* end of any combination *\/ */ } /* end of nres */ /* hBijx(p, bage, fage); */ /* fclose(ficrespijb); */ @@ -11608,7 +11763,7 @@ int back_prevalence_limit(double *p, dou int hPijx(double *p, int bage, int fage){ /*------------- h Pij x at various ages ------------*/ - + /* to be optimized with precov */ int stepsize; int agelim; int hstepm; @@ -11618,74 +11773,77 @@ int hPijx(double *p, int bage, int fage) double agedeb; double ***p3mat; - strcpy(filerespij,"PIJ_"); strcat(filerespij,fileresu); - if((ficrespij=fopen(filerespij,"w"))==NULL) { - printf("Problem with Pij resultfile: %s\n", filerespij); return 1; - fprintf(ficlog,"Problem with Pij resultfile: %s\n", filerespij); return 1; - } - printf("Computing pij: result on file '%s' \n", filerespij); - fprintf(ficlog,"Computing pij: result on file '%s' \n", filerespij); - - stepsize=(int) (stepm+YEARM-1)/YEARM; - /*if (stepm<=24) stepsize=2;*/ - - agelim=AGESUP; - hstepm=stepsize*YEARM; /* Every year of age */ - hstepm=hstepm/stepm; /* Typically 2 years, = 2/6 months = 4 */ - - /* hstepm=1; aff par mois*/ - pstamp(ficrespij); - fprintf(ficrespij,"#****** h Pij x Probability to be in state j at age x+h being in i at x "); - i1= pow(2,cptcoveff); - /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */ - /* /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */ - /* k=k+1; */ - for(nres=1; nres <= nresult; nres++) /* For each resultline */ - for(k=1; k<=i1;k++){ - if(i1 != 1 && TKresult[nres]!= k) - continue; - fprintf(ficrespij,"\n#****** "); - for(j=1;j<=cptcoveff;j++) - fprintf(ficrespij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */ - printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); - fprintf(ficrespij," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); - } - fprintf(ficrespij,"******\n"); - - for (agedeb=fage; agedeb>=bage; agedeb--){ /* If stepm=6 months */ - nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm); /* Typically 20 years = 20*12/6=40 */ - nhstepm = nhstepm/hstepm; /* Typically 40/4=10 */ - - /* nhstepm=nhstepm*YEARM; aff par mois*/ - - p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); - oldm=oldms;savm=savms; - hpxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k, nres); - fprintf(ficrespij,"# Cov Agex agex+h hpijx with i,j="); + strcpy(filerespij,"PIJ_"); strcat(filerespij,fileresu); + if((ficrespij=fopen(filerespij,"w"))==NULL) { + printf("Problem with Pij resultfile: %s\n", filerespij); return 1; + fprintf(ficlog,"Problem with Pij resultfile: %s\n", filerespij); return 1; + } + printf("Computing pij: result on file '%s' \n", filerespij); + fprintf(ficlog,"Computing pij: result on file '%s' \n", filerespij); + + stepsize=(int) (stepm+YEARM-1)/YEARM; + /*if (stepm<=24) stepsize=2;*/ + + agelim=AGESUP; + hstepm=stepsize*YEARM; /* Every year of age */ + hstepm=hstepm/stepm; /* Typically 2 years, = 2/6 months = 4 */ + + /* hstepm=1; aff par mois*/ + pstamp(ficrespij); + fprintf(ficrespij,"#****** h Pij x Probability to be in state j at age x+h being in i at x "); + i1= pow(2,cptcoveff); + /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */ + /* /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */ + /* k=k+1; */ + for(nres=1; nres <= nresult; nres++){ /* For each resultline */ + k=TKresult[nres]; + if(TKresult[nres]==0) k=1; /* To be checked for noresult */ + /* for(k=1; k<=i1;k++){ */ + /* if(i1 != 1 && TKresult[nres]!= k) */ + /* continue; */ + fprintf(ficrespij,"\n#****** "); + for(j=1;j<=cptcovs;j++){ + fprintf(ficrespij," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + /* fprintf(ficrespij,"@wV%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* for (k4=1; k4<= nsq; k4++){ /\* For each selected (single) quantitative value *\/ */ + /* printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); */ + /* fprintf(ficrespij," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); */ + } + fprintf(ficrespij,"******\n"); + + for (agedeb=fage; agedeb>=bage; agedeb--){ /* If stepm=6 months */ + nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm); /* Typically 20 years = 20*12/6=40 */ + nhstepm = nhstepm/hstepm; /* Typically 40/4=10 */ + + /* nhstepm=nhstepm*YEARM; aff par mois*/ + + p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); + oldm=oldms;savm=savms; + hpxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k, nres); + fprintf(ficrespij,"# Cov Agex agex+h hpijx with i,j="); + for(i=1; i<=nlstate;i++) + for(j=1; j<=nlstate+ndeath;j++) + fprintf(ficrespij," %1d-%1d",i,j); + fprintf(ficrespij,"\n"); + for (h=0; h<=nhstepm; h++){ + /*agedebphstep = agedeb + h*hstepm/YEARM*stepm;*/ + fprintf(ficrespij,"%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(ficrespij," %1d-%1d",i,j); - fprintf(ficrespij,"\n"); - for (h=0; h<=nhstepm; h++){ - /*agedebphstep = agedeb + h*hstepm/YEARM*stepm;*/ - fprintf(ficrespij,"%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(ficrespij," %.5f", p3mat[i][j][h]); - fprintf(ficrespij,"\n"); - } - free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); + fprintf(ficrespij," %.5f", p3mat[i][j][h]); fprintf(ficrespij,"\n"); } - /*}*/ + free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); + fprintf(ficrespij,"\n"); } - return 0; + } + /*}*/ + return 0; } int hBijx(double *p, int bage, int fage, double ***prevacurrent){ /*------------- h Bij x at various ages ------------*/ - + /* To be optimized with precov */ int stepsize; /* int agelim; */ int ageminl; @@ -11720,54 +11878,58 @@ int hPijx(double *p, int bage, int fage) /* /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */ /* k=k+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(i1 != 1 && 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,TnsdVar[Tvaraff[j]])]); - for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */ - fprintf(ficrespijb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - } - fprintf(ficrespijb,"******\n"); - if(invalidvarcomb[k]){ /* Is it necessary here? */ - 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+0.1)-1; /* Typically 20 years = 20*12/6=40 or 55*12/24=27.5-1.1=>27 */ - nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 or 28*/ - - /* nhstepm=nhstepm*YEARM; aff par mois*/ - - p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); /* We can't have it at an upper level because of nhstepm */ - /* and memory limitations if stepm is small */ - - /* 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, nres);/* Bug valgrind */ - /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm,oldm,savm, dnewm, doldm, dsavm, k); */ - fprintf(ficrespijb,"# Cov Agex agex-h hbijx with i,j="); + k=TKresult[nres]; + if(TKresult[nres]==0) k=1; /* To be checked for noresult */ + /* for(k=1; k<=i1;k++){ /\* For any combination of dummy covariates, fixed and varying *\/ */ + /* if(i1 != 1 && TKresult[nres]!= k) */ + /* continue; */ + fprintf(ficrespijb,"\n#****** "); + for(j=1;j<=cptcovs;j++){ + fprintf(ficrespijb," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]); + /* for(j=1;j<=cptcoveff;j++) */ + /* fprintf(ficrespijb,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ + /* for (j=1; j<= nsq; j++){ /\* For each selected (single) quantitative value *\/ */ + /* fprintf(ficrespijb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + } + fprintf(ficrespijb,"******\n"); + if(invalidvarcomb[k]){ /* Is it necessary here? */ + 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+0.1)-1; /* Typically 20 years = 20*12/6=40 or 55*12/24=27.5-1.1=>27 */ + nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 or 28*/ + + /* nhstepm=nhstepm*YEARM; aff par mois*/ + + p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); /* We can't have it at an upper level because of nhstepm */ + /* and memory limitations if stepm is small */ + + /* 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, nres);/* Bug valgrind */ + /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm,oldm,savm, dnewm, doldm, dsavm, k); */ + fprintf(ficrespijb,"# Cov Agex agex-h hbijx 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(i=1; i<=nlstate;i++) for(j=1; j<=nlstate+ndeath;j++) - fprintf(ficrespijb," %1d-%1d",i,j); + fprintf(ficrespijb," %.5f", p3mat[i][j][h]);/* Bug valgrind */ 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]);/* Bug valgrind */ - fprintf(ficrespijb,"\n"); - } - free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); - fprintf(ficrespijb,"\n"); - } /* end age deb */ - } /* end combination */ + } + free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm); + fprintf(ficrespijb,"\n"); + } /* end age deb */ + /* } /\* end combination *\/ */ } /* end nres */ return 0; } /* hBijx */ @@ -12163,7 +12325,7 @@ int main(int argc, char *argv[]) strcpy(model,modeltemp); } } - /* printf(" model=1+age%s modeltemp= %s, model=%s\n",model, modeltemp, model);fflush(stdout); */ + /* printf(" model=1+age%s modeltemp= %s, model=1+age+%s\n",model, modeltemp, model);fflush(stdout); */ printf("model=1+age+%s\n",model);fflush(stdout); fprintf(ficparo,"model=1+age+%s\n",model);fflush(stdout); fprintf(ficres,"model=1+age+%s\n",model);fflush(stdout); @@ -12447,7 +12609,7 @@ Please run with mle=-1 to get a correct mint=matrix(1,maxwav,firstobs,lastobs); anint=matrix(1,maxwav,firstobs,lastobs); s=imatrix(1,maxwav+1,firstobs,lastobs); /* s[i][j] health state for wave i and individual j */ - printf("BUG ncovmodel=%d NCOVMAX=%d 2**ncovmodel=%f BUG\n",ncovmodel,NCOVMAX,pow(2,ncovmodel)); + /* printf("BUG ncovmodel=%d NCOVMAX=%d 2**ncovmodel=%f BUG\n",ncovmodel,NCOVMAX,pow(2,ncovmodel)); */ tab=ivector(1,NCOVMAX); ncodemax=ivector(1,NCOVMAX); /* Number of code per covariate; if O and 1 only, 2**ncov; V1+V2+V3+V4=>16 */ ncodemaxwundef=ivector(1,NCOVMAX); /* Number of code per covariate; if - 1 O and 1 only, 2**ncov; V1+V2+V3+V4=>16 */ @@ -12594,8 +12756,8 @@ Please run with mle=-1 to get a correct } ncovcombmax=pow(2,cptcoveff); - invalidvarcomb=ivector(1, ncovcombmax); - for(i=1;iDatafile=%s Firstpass=%d La fprintf(fichtm,"
\n\ IMaCh-%s
%s \
\n\ -This file: %sTitle=%s
Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s
\n\ +This file: %s
Title=%s
Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s
\n\ \n\
\

  • Parameter files

    \n\ @@ -13742,8 +13904,8 @@ Please run with mle=-1 to get a correct printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); } for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */ - printf(" V%d=%f ",TvarsQ[j], TinvDoQresult[nres][TvarsQ[j]]); /* TvarsQ[j] gives the name of the jth quantitative (fixed or time v) */ - fprintf(ficreseij,"V%d=%f ",TvarsQ[j], TinvDoQresult[nres][TvarsQ[j]]); + printf(" V%d=%lg ",TvarsQ[j], TinvDoQresult[nres][TvarsQ[j]]); /* TvarsQ[j] gives the name of the jth quantitative (fixed or time v) */ + fprintf(ficreseij,"V%d=%lg ",TvarsQ[j], TinvDoQresult[nres][TvarsQ[j]]); } fprintf(ficreseij,"******\n"); printf("******\n"); @@ -13761,7 +13923,7 @@ Please run with mle=-1 to get a correct /*---------- State-specific expectancies and variances ------------*/ - + /* Should be moved in a function */ strcpy(filerest,"T_"); strcat(filerest,fileresu); if((ficrest=fopen(filerest,"w"))==NULL) { @@ -13835,9 +13997,9 @@ Please run with mle=-1 to get a correct printf("\n j=%d In computing T_ Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=%d cptcovs=%d, cptcoveff=%d Fixed[modelresult[nres][j]]=%d\n", j, nres, j, Dummy[modelresult[nres][j]],nres,j,modelresult[nres][j],cptcovs, cptcoveff,Fixed[modelresult[nres][j]]); /* end if dummy or quanti */ if(Dummy[modelresult[nres][j]]==0){/* Dummy variable of the variable in position modelresult in the model corresponding to j in resultline */ - printf("V%d=%d ",Tvresult[nres][j],Tresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ - fprintf(ficlog,"V%d=%d ",Tvresult[nres][j],Tresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ - fprintf(ficrest,"V%d=%d ",Tvresult[nres][j],Tresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + printf("V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + fprintf(ficlog,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + fprintf(ficrest,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ if(Fixed[modelresult[nres][j]]==0){ /* Fixed */ printf("fixed ");fprintf(ficlog,"fixed ");fprintf(ficrest,"fixed "); }else{ @@ -13847,7 +14009,9 @@ Please run with mle=-1 to get a correct /* fprintf(ficlog,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ }else if(Dummy[modelresult[nres][j]]==1){ /* Quanti variable */ /* For each selected (single) quantitative value */ - printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]); + printf(" V%d=%lg ",Tvqresult[nres][j],Tqresult[nres][j]); + fprintf(ficlog," V%d=%lg ",Tvqresult[nres][j],Tqresult[nres][j]); + fprintf(ficrest," V%d=%lg ",Tvqresult[nres][j],Tqresult[nres][j]); if(Fixed[modelresult[nres][j]]==0){ /* Fixed */ printf("fixed ");fprintf(ficlog,"fixed ");fprintf(ficrest,"fixed "); }else{ @@ -13870,14 +14034,16 @@ Please run with mle=-1 to get a correct fprintf(ficresstdeij,"\n#****** "); fprintf(ficrescveij,"\n#****** "); + /* It could have been: for(j=1;j<=cptcoveff;j++) {printf("V=%d=%lg",Tvresult[nres][cpt],TinvDoQresult[nres][Tvresult[nres][cpt]]);} */ + /* But it won't be sorted and depends on how the resultline is ordered */ for(j=1;j<=cptcoveff;j++) { fprintf(ficresstdeij,"V%d=%d ",Tvresult[nres][j],Tresult[nres][j]); /* fprintf(ficresstdeij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ /* fprintf(ficrescveij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */ } for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value, TvarsQind gives the position of a quantitative in model equation */ - fprintf(ficresstdeij," V%d=%f ",Tvar[TvarsQind[j]],Tqresult[nres][resultmodel[nres][TvarsQind[j]]]); - fprintf(ficrescveij," V%d=%f ",Tvar[TvarsQind[j]],Tqresult[nres][resultmodel[nres][TvarsQind[j]]]); + fprintf(ficresstdeij," V%d=%lg ",Tvar[TvarsQind[j]],Tqresult[nres][resultmodel[nres][TvarsQind[j]]]); + fprintf(ficrescveij," V%d=%lg ",Tvar[TvarsQind[j]],Tqresult[nres][resultmodel[nres][TvarsQind[j]]]); } fprintf(ficresstdeij,"******\n"); fprintf(ficrescveij,"******\n"); @@ -13889,7 +14055,7 @@ Please run with mle=-1 to get a correct /* fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[TnsdVar[Tvaraff[j]]])]); */ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */ /* fprintf(ficresvij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]); /\* To solve *\/ */ - fprintf(ficresvij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); /* Solved */ + fprintf(ficresvij," V%d=%lg ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); /* Solved */ } fprintf(ficresvij,"******\n"); @@ -14046,7 +14212,7 @@ Please run with mle=-1 to get a correct free_ivector(Tposprod,1,NCOVMAX); free_ivector(Tprod,1,NCOVMAX); free_ivector(Tvaraff,1,NCOVMAX); - free_ivector(invalidvarcomb,1,ncovcombmax); + free_ivector(invalidvarcomb,0,ncovcombmax); free_ivector(Tage,1,NCOVMAX); free_ivector(Tmodelind,1,NCOVMAX); free_ivector(TmodelInvind,1,NCOVMAX);