--- imach/src/imach.c 2022/08/21 09:10:30 1.333 +++ imach/src/imach.c 2022/09/04 17:40:33 1.338 @@ -1,6 +1,29 @@ -/* $Id: imach.c,v 1.333 2022/08/21 09:10:30 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... + + Revision 1.334 2022/08/25 09:08:41 brouard + Summary: In progress for quantitative + Revision 1.333 2022/08/21 09:10:30 brouard * src/imach.c (Module): Version 0.99r33 A lot of changes in reassigning covariates: my first idea was that people will always @@ -1278,25 +1301,25 @@ typedef struct { #define ODIRSEPARATOR '\\' #endif -/* $Id: imach.c,v 1.333 2022/08/21 09:10:30 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.333 $ $Date: 2022/08/21 09:10:30 $"; +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 */ int nagesqr=0, nforce=0; /* nagesqr=1 if model is including age*age, number of forces */ /* Number of covariates model (1)=V2+V1+ V3*age+V2*V4 */ /* Model(2) V1 + V2 + V3 + V8 + V7*V8 + V5*V6 + V8*age + V3*age + age*age */ -int cptcovn=0; /**< cptcovn decodemodel: number of covariates k of the models excluding age*products =6 and age*age */ +int cptcovn=0; /**< cptcovn decodemodel: number of covariates k of the models excluding age*products =6 and age*age but including products */ int cptcovt=0; /**< cptcovt: total number of covariates of the model (2) nbocc(+)+1 = 8 excepting constant and age and age*age */ -int cptcovs=0; /**< cptcovs number of simple covariates in the model V2+V1 =2 (dummy or quantit or time varying) */ -int cptcovsnq=0; /**< cptcovsnq number of simple covariates in the model but non quantitative V2+V1 =2 */ +int cptcovs=0; /**< cptcovs number of SIMPLE covariates in the model V2+V1 =2 (dummy or quantit or time varying) */ +int cptcovsnq=0; /**< cptcovsnq number of SIMPLE covariates in the model but non quantitative V2+V1 =2 */ int cptcovage=0; /**< Number of covariates with age: V3*age only =1 */ int cptcovprodnoage=0; /**< Number of covariate products without age */ -int cptcoveff=0; /* Total number of covariates to vary for printing results (2**cptcoveff combinations of dummies)(computed in tricode as cptcov) */ +int cptcoveff=0; /* Total number of single dummy covariates (fixed or time varying) to vary for printing results (2**cptcoveff combinations of dummies)(computed in tricode as cptcov) */ int ncovf=0; /* Total number of effective fixed covariates (dummy or quantitative) in the model */ int ncovv=0; /* Total number of effective (wave) varying covariates (dummy or quantitative) in the model */ int ncova=0; /* Total number of effective (wave and stepm) varying with age covariates (dummy of quantitative) in the model */ @@ -1307,6 +1330,7 @@ int nqfveff=0; /**< nqfveff Number of Qu int ntveff=0; /**< ntveff number of effective time varying variables */ int nqtveff=0; /**< ntqveff number of effective time varying quantitative variables */ int cptcov=0; /* Working variable */ +int firstobs=1, lastobs=10; /* nobs = lastobs-firstobs+1 declared globally ;*/ int nobs=10; /* Number of observations in the data lastobs-firstobs */ int ncovcombmax=NCOVMAX; /* Maximum calculated number of covariate combination = pow(2, cptcoveff) */ int npar=NPARMAX; /* Number of parameters (nlstate+ndeath-1)*nlstate*ncovmodel; */ @@ -1498,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 */ @@ -1524,12 +1549,13 @@ int *TvarsQind; #define MAXRESULTLINESPONE 10+1 int nresult=0; int parameterline=0; /* # of the parameter (type) line */ -int TKresult[MAXRESULTLINESPONE]; -int resultmodel[MAXRESULTLINESPONE][NCOVMAX];/* resultmodel[k1]=k3: k1th position in the model correspond to the k3 position in the resultline */ -int Tresult[MAXRESULTLINESPONE][NCOVMAX];/* For dummy variable , value (output) */ +int TKresult[MAXRESULTLINESPONE]; /* TKresult[nres]=k for each resultline nres give the corresponding combination of dummies */ +int resultmodel[MAXRESULTLINESPONE][NCOVMAX];/* resultmodel[k1]=k3: k1th position in the model corresponds to the k3 position in the resultline */ +int modelresult[MAXRESULTLINESPONE][NCOVMAX];/* modelresult[k3]=k1: k1th position in the model corresponds to the k3 position in the resultline */ +int Tresult[MAXRESULTLINESPONE][NCOVMAX];/* Tresult[nres][result_position]= value of the dummy variable at the result_position in the nres resultline */ int Tinvresult[MAXRESULTLINESPONE][NCOVMAX];/* Tinvresult[nres][Name of a dummy variable]= value of the variable in the result line */ double TinvDoQresult[MAXRESULTLINESPONE][NCOVMAX];/* TinvDoQresult[nres][Name of a Dummy or Q variable]= value of the variable in the result line */ -int Tvresult[MAXRESULTLINESPONE][NCOVMAX]; /* For dummy variable , variable # (output) */ +int Tvresult[MAXRESULTLINESPONE][NCOVMAX]; /* Tvresult[nres][result_position]= name of the dummy variable at the result_position in the nres resultline */ double Tqresult[MAXRESULTLINESPONE][NCOVMAX]; /* Tqresult[nres][result_position]= value of the variable at the result_position in the nres resultline */ double Tqinvresult[MAXRESULTLINESPONE][NCOVMAX]; /* For quantitative variable , value (output) */ int Tvqresult[MAXRESULTLINESPONE][NCOVMAX]; /* Tvqresult[nres][result_position]= id of the variable at the result_position in the nres resultline */ @@ -2491,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; @@ -2499,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++) { @@ -2803,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 @@ -3500,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. @@ -3831,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]);*/ @@ -3863,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] @@ -3883,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]; @@ -3893,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 */ @@ -4167,7 +4202,7 @@ double func( double *x) double funcone( double *x) { /* Same as func but slower because of a lot of printf and if */ - int i, ii, j, k, mi, d, kk; + int i, ii, j, k, mi, d, kk, kf=0; int ioffset=0; double l, ll[NLSTATEMAX+1], cov[NCOVMAX+1]; double **out; @@ -4190,13 +4225,18 @@ 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 */ /* for (k=1; k<=cptcovn;k++) cov[2+nagesqr+k]=covar[Tvar[k]][i]; */ /* for (k=1; k<=ncoveff;k++){ /\* Simple and product fixed Dummy covariates without age* products *\/ */ - for (k=1; k<=ncovf;k++){ /* Simple and product fixed covariates without age* products *//* Missing values are set to -1 but should be dropped */ - cov[ioffset+TvarFind[k]]=covar[Tvar[TvarFind[k]]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V1 is fixed (k=6)*/ + for (kf=1; kf<=ncovf;kf++){ /* Simple and product fixed covariates without age* products *//* Missing values are set to -1 but should be dropped */ + cov[ioffset+TvarFind[kf]]=covar[Tvar[TvarFind[kf]]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V1 is fixed (k=6)*/ /* cov[ioffset+TvarFind[1]]=covar[Tvar[TvarFind[1]]][i]; */ /* cov[2+6]=covar[Tvar[6]][i]; */ /* cov[2+6]=covar[2][i]; V2 */ @@ -4207,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?)*\/ */ /* } */ @@ -4264,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){ */ @@ -4296,27 +4361,33 @@ double funcone( double *x) ipmx +=1; sw += weight[i]; ll[s[mw[mi][i]][i]] += 2*weight[i]*lli; - /*printf("i=%6d s1=%1d s2=%1d mi=%1d mw=%1d dh=%3d prob=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],out[s1][s2],savm[s1][s2]); */ + /* printf("Funcone i=%6d s1=%1d s2=%1d mi=%1d mw=%1d dh=%3d prob=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); */ if(globpr){ fprintf(ficresilk,"%09ld %6.1f %6.1f %6d %2d %2d %2d %2d %3d %15.6f %8.4f %8.3f\ %11.6f %11.6f %11.6f ", \ num[i], agebegin, ageend, i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],weight[i]*gipmx/gsw, 2*weight[i]*lli,(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); + /* printf("%09ld %6.1f %6.1f %6d %2d %2d %2d %2d %3d %15.6f %8.4f %8.3f\ */ + /* %11.6f %11.6f %11.6f ", \ */ + /* num[i], agebegin, ageend, i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],weight[i]*gipmx/gsw, */ + /* 2*weight[i]*lli,(s2==-1? -1: out[s1][s2]),(s2==-1? -1: savm[s1][s2])); */ for(k=1,llt=0.,l=0.; k<=nlstate; k++){ llt +=ll[k]*gipmx/gsw; fprintf(ficresilk," %10.6f",-ll[k]*gipmx/gsw); + /* printf(" %10.6f",-ll[k]*gipmx/gsw); */ } fprintf(ficresilk," %10.6f\n", -llt); + /* printf(" %10.6f\n", -llt); */ } - } /* end of wave */ -} /* end of individual */ -for(k=1,l=0.; k<=nlstate; k++) l += ll[k]; + } /* end of wave */ + } /* end of individual */ + for(k=1,l=0.; k<=nlstate; k++) l += ll[k]; /* printf("l1=%f l2=%f ",ll[1],ll[2]); */ -l= l*ipmx/sw; /* To get the same order of magnitude as if weight=1 for every body */ -if(globpr==0){ /* First time we count the contributions and weights */ - gipmx=ipmx; - gsw=sw; -} + l= l*ipmx/sw; /* To get the same order of magnitude as if weight=1 for every body */ + if(globpr==0){ /* First time we count the contributions and weights */ + gipmx=ipmx; + gsw=sw; + } return -l; } @@ -4975,7 +5046,7 @@ Title=%s
Datafile=%s Firstpass=%d La j1=0; /* j=ncoveff; /\* Only fixed dummy covariates *\/ */ - j=cptcoveff; /* Only dummy covariates of the model */ + j=cptcoveff; /* Only simple dummy covariates used in the model */ /* j=cptcovn; /\* Only dummy covariates of the model *\/ */ if (cptcovn<1) {j=1;ncodemax[1]=1;} @@ -4996,7 +5067,7 @@ Title=%s
Datafile=%s Firstpass=%d La /* if a constant only model, one pass to compute frequency tables and to write it on ficresp */ /* Loop on nj=1 or 2 if dummy covariates j!=0 - * Loop on j1(1 to 2**cptcovn) covariate combination + * Loop on j1(1 to 2**cptcoveff) covariate combination * freq[s1][s2][iage] =0. * Loop on iind * ++freq[s1][s2][iage] weighted @@ -5021,7 +5092,7 @@ Title=%s
Datafile=%s Firstpass=%d La if(nj==1) j=0; /* First pass for the constant */ else{ - j=cptcovs; /* Other passes for the covariate values */ + j=cptcoveff; /* Other passes for the covariate values number of simple covariates in the model V2+V1 =2 (simple dummy fixed or time varying) */ } first=1; for (j1 = 1; j1 <= (int) pow(2,j); j1++){ /* Loop on all dummy covariates combination of the model, ie excluding quantitatives, V4=0, V3=0 for example, fixed or varying covariates */ @@ -5057,13 +5128,16 @@ Title=%s
Datafile=%s Firstpass=%d La bool=1; if(j !=0){ if(anyvaryingduminmodel==0){ /* If All fixed covariates */ - if (cptcovn >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */ - for (z1=1; z1<=cptcovn; z1++) { /* loops on covariates in the model */ + if (cptcoveff >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */ + for (z1=1; z1<=cptcoveff; z1++) { /* loops on covariates in the model */ /* if(Tvaraff[z1] ==-20){ */ /* /\* sumnew+=cotvar[mw[mi][iind]][z1][iind]; *\/ */ /* }else if(Tvaraff[z1] ==-10){ */ /* /\* sumnew+=coqvar[z1][iind]; *\/ */ /* }else */ /* TODO TODO codtabm(j1,z1) or codtabm(j1,Tvaraff[z1]]z1)*/ + /* if( iind >=imx-3) printf("Searching error iind=%d Tvaraff[z1]=%d covar[Tvaraff[z1]][iind]=%.f TnsdVar[Tvaraff[z1]]=%d, cptcoveff=%d, cptcovs=%d \n",iind, Tvaraff[z1], covar[Tvaraff[z1]][iind],TnsdVar[Tvaraff[z1]],cptcoveff, cptcovs); */ + if(Tvaraff[z1]<1 || Tvaraff[z1]>=NCOVMAX) + printf("Error Tvaraff[z1]=%d<1 or >=%d, cptcoveff=%d model=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 */ @@ -5073,7 +5147,7 @@ Title=%s
Datafile=%s Firstpass=%d La /* For j1=7 in V1+V2+V3+V4 = 0 1 1 0 and codtabm(7,3)=1 and nbcde[3][?]=1*/ } /* Onlyf fixed */ } /* end z1 */ - } /* cptcovn > 0 */ + } /* cptcoveff > 0 */ } /* end any */ }/* end j==0 */ if (bool==1){ /* We selected an individual iind satisfying combination j1 (V4=1 V3=0) or all fixed covariates */ @@ -5082,7 +5156,7 @@ Title=%s
Datafile=%s Firstpass=%d La m=mw[mi][iind]; if(j!=0){ if(anyvaryingduminmodel==1){ /* Some are varying covariates */ - for (z1=1; z1<=cptcovn; z1++) { + for (z1=1; z1<=cptcoveff; z1++) { if( Fixed[Tmodelind[z1]]==1){ iv= Tvar[Tmodelind[z1]]-ncovcol-nqv; if (cotvar[m][iv][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]) /* iv=1 to ntv, right modality. If covariate's @@ -5090,7 +5164,12 @@ Title=%s
Datafile=%s Firstpass=%d La constant and age model which counts them. */ bool=0; /* not selected */ }else if( Fixed[Tmodelind[z1]]== 0) { /* fixed */ - if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]) { + /* i1=Tvaraff[z1]; */ + /* i2=TnsdVar[i1]; */ + /* i3=nbcode[i1][i2]; */ + /* i4=covar[i1][iind]; */ + /* if(i4 != i3){ */ + if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]) { /* Bug valgrind */ bool=0; } } @@ -5151,9 +5230,9 @@ Title=%s
Datafile=%s Firstpass=%d La /* fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/ - if(cptcovn==0 && nj==1) /* no covariate and first pass */ + if(cptcoveff==0 && nj==1) /* no covariate and first pass */ pstamp(ficresp); - if (cptcovn>0 && j!=0){ + if (cptcoveff>0 && j!=0){ pstamp(ficresp); printf( "\n#********** Variable "); fprintf(ficresp, "\n#********** Variable "); @@ -5206,14 +5285,17 @@ Title=%s
Datafile=%s Firstpass=%d La /* } */ fprintf(ficresphtm,""); - if((cptcovn==0 && nj==1)|| nj==2 ) /* no covariate and first pass */ + if((cptcoveff==0 && nj==1)|| nj==2 ) /* no covariate and first pass */ fprintf(ficresp, " Age"); - if(nj==2) for (z1=1; z1<=cptcovn; z1++) fprintf(ficresp, " V%d=%d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); + if(nj==2) for (z1=1; z1<=cptcoveff; z1++) { + printf(" V%d=%d, z1=%d, Tvaraff[z1]=%d, j1=%d, TnsdVar[Tvaraff[%d]]=%d |",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])], z1, Tvaraff[z1], j1,z1,TnsdVar[Tvaraff[z1]]); + fprintf(ficresp, " V%d=%d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); + } for(i=1; i<=nlstate;i++) { - if((cptcovn==0 && nj==1)|| nj==2 ) fprintf(ficresp," Prev(%d) N(%d) N ",i,i); + if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp," Prev(%d) N(%d) N ",i,i); fprintf(ficresphtm, "",i,i); } - if((cptcovn==0 && nj==1)|| nj==2 ) fprintf(ficresp, "\n"); + if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp, "\n"); fprintf(ficresphtm, "\n"); /* Header of frequency table by age */ @@ -5281,14 +5363,14 @@ Title=%s
Datafile=%s Firstpass=%d La } /* Writing ficresp */ - if(cptcovn==0 && nj==1){ /* no covariate and first pass */ + if(cptcoveff==0 && nj==1){ /* no covariate and first pass */ if( iage <= iagemax){ fprintf(ficresp," %d",iage); } }else if( nj==2){ if( iage <= iagemax){ fprintf(ficresp," %d",iage); - for (z1=1; z1<=cptcovn; z1++) fprintf(ficresp, " %d %d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); + for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresp, " %d %d",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); } } for(s1=1; s1 <=nlstate ; s1++){ @@ -5303,7 +5385,7 @@ Title=%s
Datafile=%s Firstpass=%d La } if( iage <= iagemax){ if(pos>=1.e-5){ - if(cptcovn==0 && nj==1){ /* no covariate and first pass */ + if(cptcoveff==0 && nj==1){ /* no covariate and first pass */ fprintf(ficresp," %.5f %.0f %.0f",prop[s1][iage]/pospropta, prop[s1][iage],pospropta); }else if( nj==2){ fprintf(ficresp," %.5f %.0f %.0f",prop[s1][iage]/pospropta, prop[s1][iage],pospropta); @@ -5312,7 +5394,7 @@ Title=%s
Datafile=%s Firstpass=%d La /*probs[iage][s1][j1]= pp[s1]/pos;*/ /*printf("\niage=%d s1=%d j1=%d %.5f %.0f %.0f %f",iage,s1,j1,pp[s1]/pos, pp[s1],pos,probs[iage][s1][j1]);*/ } else{ - if((cptcovn==0 && nj==1)|| nj==2 ) fprintf(ficresp," NaNq %.0f %.0f",prop[s1][iage],pospropta); + if((cptcoveff==0 && nj==1)|| nj==2 ) fprintf(ficresp," NaNq %.0f %.0f",prop[s1][iage],pospropta); fprintf(ficresphtm,"",iage, prop[s1][iage],pospropta); } } @@ -5338,7 +5420,7 @@ Title=%s
Datafile=%s Firstpass=%d La } fprintf(ficresphtmfr,"\n "); fprintf(ficresphtm,"\n"); - if((cptcovn==0 && nj==1)|| nj==2 ) { + if((cptcoveff==0 && nj==1)|| nj==2 ) { if(iage <= iagemax) fprintf(ficresp,"\n"); } @@ -5365,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,"
AgePrev(%d)N(%d)N%dNaNq%.0f%.0f
\n"); @@ -5596,7 +5678,7 @@ void prevalence(double ***probs, double if (cptcovn<1) {j=1;ncodemax[1]=1;} first=0; - for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of covariate */ + for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of simple dummy covariates */ for (i=1; i<=nlstate; i++) for(iage=iagemin-AGEMARGE; iage <= iagemax+4+AGEMARGE; iage++) prop[i][iage]=0.0; @@ -5920,7 +6002,7 @@ void concatwav(int wav[], int **dh, int nbcode[k][j]=0; /* Valgrind */ /* Loop on covariates without age and products and no quantitative variable */ - for (k=1; k<=cptcovt; k++) { /* From model V1 + V2*age + V3 + V3*V4 keeps V1 + V3 = 2 only */ + for (k=1; k<=cptcovt; k++) { /* cptcovt: total number of covariates of the model (2) nbocc(+)+1 = 8 excepting constant and age and age*age */ for (j=-1; (j < maxncov); j++) Ndum[j]=0; if(Dummy[k]==0 && Typevar[k] !=1){ /* Dummy covariate and not age product */ switch(Fixed[k]) { @@ -6018,8 +6100,12 @@ void concatwav(int wav[], int **dh, int break; } /* end switch */ } /* end dummy test */ - if(Dummy[k]==1 && Typevar[k] !=1){ /* Dummy covariate and not age product */ + if(Dummy[k]==1 && Typevar[k] !=1){ /* Quantitative covariate and not age product */ for (i=1; i<=imx; i++) { /* Loop on individuals: reads the data file to get the maximum value of the modality of this covariate Vj*/ + if(Tvar[k]<=0 || Tvar[k]>=NCOVMAX){ + printf("Error k=%d \n",k); + exit(1); + } if(isnan(covar[Tvar[k]][i])){ printf("ERROR, IMaCh doesn't treat fixed quantitative covariate with missing values V%d=., individual %d will be skipped.\n",Tvar[k],i); fprintf(ficlog,"ERROR, currently IMaCh doesn't treat covariate with missing values V%d=., individual %d will be skipped.\n",Tvar[k],i); @@ -6027,7 +6113,7 @@ void concatwav(int wav[], int **dh, int exit(1); } } - } + } /* end Quanti */ } /* end of loop on model-covariate k. nbcode[Tvark][1]=-1, nbcode[Tvark][1]=0 and nbcode[Tvark][2]=1 sets the value of covariate k*/ for (k=-1; k< maxncov; k++) Ndum[k]=0; @@ -6041,13 +6127,22 @@ void concatwav(int wav[], int **dh, int ij=0; /* for (i=0; i<= maxncov-1; i++) { /\* modmaxcovj is unknown here. Only Ndum[2(V2),3(age*V3), 5(V3*V2) 6(V1*V4) *\/ */ - for (k=1; k<= cptcovt; k++) { /* modmaxcovj is unknown here. Only Ndum[2(V2),3(age*V3), 5(V3*V2) 6(V1*V4) */ + for (k=1; k<= cptcovt; k++) { /* cptcovt: total number of covariates of the model (2) nbocc(+)+1 = 8 excepting constant and age and age*age */ + /* modmaxcovj is unknown here. Only Ndum[2(V2),3(age*V3), 5(V3*V2) 6(V1*V4) */ /*printf("Ndum[%d]=%d\n",i, Ndum[i]);*/ /* if((Ndum[i]!=0) && (i<=ncovcol)){ /\* Tvar[i] <= ncovmodel ? *\/ */ - if(Ndum[Tvar[k]]!=0 && Dummy[k] == 0 && Typevar[k]==0){ /* Only Dummy and non empty in the model */ + if(Ndum[Tvar[k]]!=0 && Dummy[k] == 0 && Typevar[k]==0){ /* Only Dummy simple and non empty in the model */ + /* Typevar[k] =0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product */ + /* Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product*/ /* If product not in single variable we don't print results */ /*printf("diff Ndum[%d]=%d\n",i, Ndum[i]);*/ - ++ij;/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, */ + ++ij;/* V5 + V4 + V3 + V4*V3 + V5*age + V2 + V1*V2 + V1*age + V1, *//* V5 quanti, V2 quanti, V4, V3, V1 dummies */ + /* k= 1 2 3 4 5 6 7 8 9 */ + /* Tvar[k]= 5 4 3 6 5 2 7 1 1 */ + /* ij 1 2 3 */ + /* Tvaraff[ij]= 4 3 1 */ + /* Tmodelind[ij]=2 3 9 */ + /* TmodelInvind[ij]=2 1 1 */ Tvaraff[ij]=Tvar[k]; /* For printing combination *//* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, Tvar {5, 4, 3, 6, 5, 2, 7, 1, 1} Tvaraff={4, 3, 1} V4, V3, V1*/ Tmodelind[ij]=k; /* Tmodelind: index in model of dummies Tmodelind[1]=2 V4: pos=2; V3: pos=3, V1=9 {2, 3, 9, ?, ?,} */ TmodelInvind[ij]=Tvar[k]- ncovcol-nqv; /* Inverse TmodelInvind[2=V4]=2 second dummy varying cov (V4)4-1-1 {0, 2, 1, } TmodelInvind[3]=1 */ @@ -6063,7 +6158,7 @@ void concatwav(int wav[], int **dh, int } /* Tvaraff[1]@5 {3, 4, -20, 0, 0} Very strange */ /* ij--; */ /* cptcoveff=ij; /\*Number of total covariates*\/ */ - *cptcov=ij; /* cptcov= Number of total real effective covariates: effective (used as cptcoveff in other functions) + *cptcov=ij; /* cptcov= Number of total real effective simple dummies (fixed or time arying) effective (used as cptcoveff in other functions) * because they can be excluded from the model and real * if in the model but excluded because missing values, but how to get k from ij?*/ for(j=ij+1; j<= cptcovt; j++){ @@ -6196,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; @@ -6464,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 */ - 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); @@ -7099,7 +7200,7 @@ To be simple, these graphs help to under for(nres=1;nres <=nresult; nres++){ /* For each resultline */ for(j1=1; j1<=tj;j1++){ /* For any combination of dummy covariates, fixed and varying */ - printf("Varprob TKresult[nres]=%d j1=%d, nres=%d, cptcovn=%d, cptcoveff=%d tj=%d \n", TKresult[nres], j1, nres, cptcovn, cptcoveff, tj); + printf("Varprob TKresult[nres]=%d j1=%d, nres=%d, cptcovn=%d, cptcoveff=%d tj=%d cptcovs=%d\n", TKresult[nres], j1, nres, cptcovn, cptcoveff, tj, cptcovs); if(tj != 1 && TKresult[nres]!= j1) continue; @@ -7107,25 +7208,63 @@ To be simple, these graphs help to under /* for(nres=1;nres <=1; nres++){ /\* For each resultline *\/ */ /* /\* for(nres=1;nres <=nresult; nres++){ /\\* For each resultline *\\/ *\/ */ if (cptcovn>0) { - fprintf(ficresprob, "\n#********** Variable "); - for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprob, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); - fprintf(ficresprob, "**********\n#\n"); + fprintf(ficresprob, "\n#********** Variable "); fprintf(ficresprobcov, "\n#********** Variable "); - for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprobcov, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); + fprintf(ficgp, "\n#********** Variable "); + fprintf(fichtmcov, "\n
********** Variable "); + fprintf(ficresprobcor, "\n#********** Variable "); + + /* Including quantitative variables of the resultline to be done */ + for (z1=1; z1<=cptcovs; z1++){ /* Loop on each variable of this resultline */ + printf("Varprob modelresult[%d][%d]=%d model=1+age+%s \n",nres, z1, modelresult[nres][z1], model); + 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 */ + fprintf(ficresprobcov,"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 */ + fprintf(ficgp,"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 */ + fprintf(fichtmcov,"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 */ + fprintf(ficresprobcor,"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 */ + fprintf(ficresprob,"fixed "); + fprintf(ficresprobcov,"fixed "); + fprintf(ficgp,"fixed "); + fprintf(fichtmcov,"fixed "); + fprintf(ficresprobcor,"fixed "); + }else{ + fprintf(ficresprob,"varyi "); + fprintf(ficresprobcov,"varyi "); + fprintf(ficgp,"varyi "); + fprintf(fichtmcov,"varyi "); + fprintf(ficresprobcor,"varyi "); + } + }else if(Dummy[modelresult[nres][z1]]==1){ /* Quanti variable */ + /* For each selected (single) quantitative value */ + fprintf(ficresprob," V%d=%lg ",Tvqresult[nres][z1],Tqresult[nres][z1]); + if(Fixed[modelresult[nres][z1]]==0){ /* Fixed */ + fprintf(ficresprob,"fixed "); + fprintf(ficresprobcov,"fixed "); + fprintf(ficgp,"fixed "); + fprintf(fichtmcov,"fixed "); + fprintf(ficresprobcor,"fixed "); + }else{ + fprintf(ficresprob,"varyi "); + fprintf(ficresprobcov,"varyi "); + fprintf(ficgp,"varyi "); + fprintf(fichtmcov,"varyi "); + fprintf(ficresprobcor,"varyi "); + } + }else{ + printf("Error in varprob() Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=V%d cptcovs=%d, cptcoveff=%d \n", nres, z1, Dummy[modelresult[nres][z1]],nres,z1,modelresult[nres][z1],cptcovs, cptcoveff); /* end if dummy or quanti */ + fprintf(ficlog,"Error in varprob() Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=V%d cptcovs=%d, cptcoveff=%d \n", nres, z1, Dummy[modelresult[nres][z1]],nres,z1,modelresult[nres][z1],cptcovs, cptcoveff); /* end if dummy or quanti */ + exit(1); + } + } /* End loop on variable of this resultline */ + /* for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprob, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); */ + fprintf(ficresprob, "**********\n#\n"); fprintf(ficresprobcov, "**********\n#\n"); - - fprintf(ficgp, "\n#********** Variable "); - for (z1=1; z1<=cptcoveff; z1++) fprintf(ficgp, " V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); fprintf(ficgp, "**********\n#\n"); - - - fprintf(fichtmcov, "\n
********** Variable "); - /* for (z1=1; z1<=cptcoveff; z1++) fprintf(fichtm, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); */ - for (z1=1; z1<=cptcoveff; z1++) fprintf(fichtmcov, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); fprintf(fichtmcov, "**********\n
"); - - fprintf(ficresprobcor, "\n#********** Variable "); - for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprobcor, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,TnsdVar[Tvaraff[z1]])]); fprintf(ficresprobcor, "**********\n#"); if(invalidvarcomb[j1]){ fprintf(ficgp,"\n#Combination (%d) ignored because no cases \n",j1); @@ -7137,57 +7276,66 @@ To be simple, these graphs help to under trgradg=matrix(1,(nlstate)*(nlstate+ndeath),1,npar); gp=vector(1,(nlstate)*(nlstate+ndeath)); gm=vector(1,(nlstate)*(nlstate+ndeath)); - for (age=bage; age<=fage; age ++){ + for (age=bage; age<=fage; age ++){ /* Fo each age we feed the model equation with covariates, using precov as in hpxij() ? */ cov[2]=age; if(nagesqr==1) cov[3]= age*age; - /* for (k=1; k<=cptcovn;k++) { */ - /* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,k)]; */ - for (k=1; k<=nsd;k++) { /* For single dummy covariates only */ - /* Here comes the value of the covariate 'j1' after renumbering k with single dummy covariates */ - cov[2+nagesqr+TvarsDind[k]]=nbcode[TvarsD[k]][codtabm(j1,TnsdVar[TvarsD[k]])]; - /*cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,Tvar[k])];*//* j1 1 2 3 4 - * 1 1 1 1 1 - * 2 2 1 1 1 - * 3 1 2 1 1 - */ - /* nbcode[1][1]=0 nbcode[1][2]=1;*/ - } - /* V2+V1+V4+V3*age Tvar[4]=3 ; V1+V2*age Tvar[2]=2; V1+V1*age Tvar[2]=1, Tage[1]=2 */ - /* ) p nbcode[Tvar[Tage[k]]][(1 & (ij-1) >> (k-1))+1] */ - /*for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */ - for (k=1; k<=cptcovage;k++){ /* For product with age */ - if(Dummy[Tage[k]]==2){ /* dummy with age */ - cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(j1,TnsdVar[Tvar[Tage[k]]])]*cov[2]; - /* cov[++k1]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; */ - } else if(Dummy[Tage[k]]==3){ /* quantitative with age */ - printf("Internal IMaCh error, don't know which value for quantitative covariate with age, Tage[k]%d, k=%d, Tvar[Tage[k]]=V%d, age=%d\n",Tage[k],k ,Tvar[Tage[k]], (int)cov[2]); - /* cov[2+nagesqr+Tage[k]]=meanq[k]/idq[k]*cov[2];/\* Using the mean of quantitative variable Tvar[Tage[k]] /\* Tqresult[nres][k]; *\/ */ - /* exit(1); */ - /* cov[++k1]=Tqresult[nres][k]; */ - } - /* cov[2+Tage[k]+nagesqr]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; */ - } - for (k=1; k<=cptcovprod;k++){/* For product without age */ - if(Dummy[Tvard[k][1]]==0){ - if(Dummy[Tvard[k][2]]==0){ - cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(j1,TnsdVar[Tvard[k][1]])] * nbcode[Tvard[k][2]][codtabm(j1,TnsdVar[Tvard[k][2]])]; - /* cov[++k1]=nbcode[Tvard[k][1]][codtabm(ij,k)] * nbcode[Tvard[k][2]][codtabm(ij,k)]; */ - }else{ /* Should we use the mean of the quantitative variables? */ - cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(j1,TnsdVar[Tvard[k][1]])] * Tqresult[nres][resultmodel[nres][k]]; - /* cov[++k1]=nbcode[Tvard[k][1]][codtabm(ij,k)] * Tqresult[nres][k]; */ - } + /* New code end of combination but for each resultline */ + for(k1=1;k1<=cptcovt;k1++){ /* loop on model equation (including products) */ + if(Typevar[k1]==1){ /* A product with age */ + cov[2+nagesqr+k1]=precov[nres][k1]*cov[2]; }else{ - if(Dummy[Tvard[k][2]]==0){ - cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][2]][codtabm(j1,TnsdVar[Tvard[k][2]])] * Tqinvresult[nres][TnsdVar[Tvard[k][1]]]; - /* cov[++k1]=nbcode[Tvard[k][2]][codtabm(ij,k)] * Tqinvresult[nres][Tvard[k][1]]; */ - }else{ - cov[2+nagesqr+Tprod[k]]=Tqinvresult[nres][TnsdVar[Tvard[k][1]]]* Tqinvresult[nres][TnsdVar[Tvard[k][2]]]; - /* cov[++k1]=Tqinvresult[nres][Tvard[k][1]]* Tqinvresult[nres][Tvard[k][2]]; */ - } + cov[2+nagesqr+k1]=precov[nres][k1]; } - /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)]; */ - } + }/* End of loop on model equation */ +/* Old code */ + /* /\* for (k=1; k<=cptcovn;k++) { *\/ */ + /* /\* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,k)]; *\/ */ + /* for (k=1; k<=nsd;k++) { /\* For single dummy covariates only *\/ */ + /* /\* Here comes the value of the covariate 'j1' after renumbering k with single dummy covariates *\/ */ + /* cov[2+nagesqr+TvarsDind[k]]=nbcode[TvarsD[k]][codtabm(j1,TnsdVar[TvarsD[k]])]; */ + /* /\*cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,Tvar[k])];*\//\* j1 1 2 3 4 */ + /* * 1 1 1 1 1 */ + /* * 2 2 1 1 1 */ + /* * 3 1 2 1 1 */ + /* *\/ */ + /* /\* nbcode[1][1]=0 nbcode[1][2]=1;*\/ */ + /* } */ + /* /\* V2+V1+V4+V3*age Tvar[4]=3 ; V1+V2*age Tvar[2]=2; V1+V1*age Tvar[2]=1, Tage[1]=2 *\/ */ + /* /\* ) p nbcode[Tvar[Tage[k]]][(1 & (ij-1) >> (k-1))+1] *\/ */ + /* /\*for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; *\/ */ + /* for (k=1; k<=cptcovage;k++){ /\* For product with age *\/ */ + /* if(Dummy[Tage[k]]==2){ /\* dummy with age *\/ */ + /* cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(j1,TnsdVar[Tvar[Tage[k]]])]*cov[2]; */ + /* /\* cov[++k1]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; *\/ */ + /* } else if(Dummy[Tage[k]]==3){ /\* quantitative with age *\/ */ + /* printf("Internal IMaCh error, don't know which value for quantitative covariate with age, Tage[k]%d, k=%d, Tvar[Tage[k]]=V%d, age=%d\n",Tage[k],k ,Tvar[Tage[k]], (int)cov[2]); */ + /* /\* cov[2+nagesqr+Tage[k]]=meanq[k]/idq[k]*cov[2];/\\* Using the mean of quantitative variable Tvar[Tage[k]] /\\* Tqresult[nres][k]; *\\/ *\/ */ + /* /\* exit(1); *\/ */ + /* /\* cov[++k1]=Tqresult[nres][k]; *\/ */ + /* } */ + /* /\* cov[2+Tage[k]+nagesqr]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; *\/ */ + /* } */ + /* for (k=1; k<=cptcovprod;k++){/\* For product without age *\/ */ + /* if(Dummy[Tvard[k][1]]==0){ */ + /* if(Dummy[Tvard[k][2]]==0){ */ + /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(j1,TnsdVar[Tvard[k][1]])] * nbcode[Tvard[k][2]][codtabm(j1,TnsdVar[Tvard[k][2]])]; */ + /* /\* cov[++k1]=nbcode[Tvard[k][1]][codtabm(ij,k)] * nbcode[Tvard[k][2]][codtabm(ij,k)]; *\/ */ + /* }else{ /\* Should we use the mean of the quantitative variables? *\/ */ + /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(j1,TnsdVar[Tvard[k][1]])] * Tqresult[nres][resultmodel[nres][k]]; */ + /* /\* cov[++k1]=nbcode[Tvard[k][1]][codtabm(ij,k)] * Tqresult[nres][k]; *\/ */ + /* } */ + /* }else{ */ + /* if(Dummy[Tvard[k][2]]==0){ */ + /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][2]][codtabm(j1,TnsdVar[Tvard[k][2]])] * Tqinvresult[nres][TnsdVar[Tvard[k][1]]]; */ + /* /\* cov[++k1]=nbcode[Tvard[k][2]][codtabm(ij,k)] * Tqinvresult[nres][Tvard[k][1]]; *\/ */ + /* }else{ */ + /* cov[2+nagesqr+Tprod[k]]=Tqinvresult[nres][TnsdVar[Tvard[k][1]]]* Tqinvresult[nres][TnsdVar[Tvard[k][2]]]; */ + /* /\* cov[++k1]=Tqinvresult[nres][Tvard[k][1]]* Tqinvresult[nres][Tvard[k][2]]; *\/ */ + /* } */ + /* } */ + /* /\* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)]; *\/ */ + /* } */ /* For each age and combination of dummy covariates we slightly move the parameters of delti in order to get the gradient*/ for(theta=1; theta <=npar; theta++){ for(i=1; i<=npar; i++) @@ -7439,29 +7587,39 @@ void printinghtml(char fileresu[], char jj1=0; fprintf(fichtm," \n"); fprintf(fichtm,"\ @@ -7623,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); @@ -7634,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); @@ -7715,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); @@ -7793,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"); */ @@ -7831,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*/ @@ -7930,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"); @@ -7998,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]){ @@ -8055,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"); @@ -8108,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"); @@ -8168,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"); @@ -8226,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"); @@ -8284,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"); @@ -8401,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"); @@ -8476,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 *\/ */ @@ -8553,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); @@ -8793,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); @@ -9362,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<= 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<=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]]); */ + /* } */ fprintf(ficresvpl,"******\n"); printf("******\n"); fprintf(ficlog,"******\n"); @@ -9419,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"); @@ -10265,7 +10488,7 @@ int decoderesult( char resultline[], int int j=0, k=0, k1=0, k2=0, k3=0, k4=0, match=0, k2q=0, k3q=0, k4q=0; char resultsav[MAXLINE]; /* int resultmodel[MAXLINE]; */ - int modelresult[MAXLINE]; + /* int modelresult[MAXLINE]; */ char stra[80], strb[80], strc[80], strd[80],stre[80]; removefirstspace(&resultline); @@ -10274,18 +10497,18 @@ int decoderesult( char resultline[], int strcpy(resultsav,resultline); printf("Decoderesult resultsav=\"%s\" resultline=\"%s\"\n", resultsav, resultline); if (strlen(resultsav) >1){ - j=nbocc(resultsav,'='); /**< j=Number of covariate values'=' */ + j=nbocc(resultsav,'='); /**< j=Number of covariate values'=' in this resultline */ } if(j == 0){ /* Resultline but no = */ TKresult[nres]=0; /* Combination for the nresult and the model */ return (0); } if( j != cptcovs ){ /* Be careful if a variable is in a product but not single */ - printf("ERROR: the number of variables in the resultline which is %d, differs from the number %d of variables used in the model line, %s.\n",j, cptcovs, model); - fprintf(ficlog,"ERROR: the number of variables in the resultline which is %d, differs from the number %d of variables used in the model line, %s.\n",j, cptcovs, model); + printf("ERROR: the number of variables in the resultline which is %d, differs from the number %d of single variables used in the model line, %s.\n",j, cptcovs, model); + fprintf(ficlog,"ERROR: the number of variables in the resultline which is %d, differs from the number %d of single variables used in the model line, %s.\n",j, cptcovs, model); /* return 1;*/ } - for(k=1; k<=j;k++){ /* Loop on any covariate of the result line */ + for(k=1; k<=j;k++){ /* Loop on any covariate of the RESULT LINE */ if(nbocc(resultsav,'=') >1){ cutl(stra,strb,resultsav,' '); /* keeps in strb after the first ' ' (stra is the rest of the resultline to be analyzed in the next loop *//* resultsav= "V4=1 V5=25.1 V3=0" stra= "V5=25.1 V3=0" strb= "V4=1" */ /* If resultsav= "V4= 1 V5=25.1 V3=0" with a blank then strb="V4=" and stra="1 V5=25.1 V3=0" */ @@ -10309,20 +10532,20 @@ int decoderesult( char resultline[], int } /* Checking for missing or useless values in comparison of current model needs */ /* Feeds resultmodel[nres][k1]=k2 for k1th product covariate with age in the model equation fed by the index k2 of the resutline*/ - for(k1=1; k1<= cptcovt ;k1++){ /* Loop on model. model line V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ + for(k1=1; k1<= cptcovt ;k1++){ /* Loop on MODEL LINE V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ if(Typevar[k1]==0){ /* Single covariate in model */ /* 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product */ match=0; for(k2=1; k2 <=j;k2++){/* Loop on resultline. In result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */ if(Tvar[k1]==Tvarsel[k2]) {/* Tvar is coming from the model, Tvarsel from the result. Tvar[1]=5 == Tvarsel[2]=5 */ - modelresult[k2]=k1;/* modelresult[2]=1 modelresult[1]=2 modelresult[3]=3 modelresult[6]=4 modelresult[9]=5 */ + modelresult[nres][k2]=k1;/* modelresult[2]=1 modelresult[1]=2 modelresult[3]=3 modelresult[6]=4 modelresult[9]=5 */ match=1; /* modelresult of k2 variable of resultline is identical to k1 variable of the model good */ break; } } 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*/ @@ -10330,7 +10553,7 @@ int decoderesult( char resultline[], int match=0; for(k2=1; k2 <=j;k2++){/* Loop on resultline. jth occurence of = signs in the result line. In result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */ if(Tvar[k1]==Tvarsel[k2]) {/* Tvar is coming from the model, Tvarsel from the result. Tvar[1]=5 == Tvarsel[2]=5 */ - modelresult[k2]=k1;/* we found a Vn=1 corrresponding to Vn*age in the model modelresult[2]=1 modelresult[1]=2 modelresult[3]=3 modelresult[6]=4 modelresult[9]=5 */ + modelresult[nres][k2]=k1;/* we found a Vn=1 corrresponding to Vn*age in the model modelresult[2]=1 modelresult[1]=2 modelresult[3]=3 modelresult[6]=4 modelresult[9]=5 */ resultmodel[nres][k1]=k2; /* Added here */ printf("Decoderesult first modelresult[k2=%d]=%d (k1) V%d*AGE\n",k2,k1,Tvar[k1]); match=1; /* modelresult of k2 variable of resultline is identical to k1 variable of the model good */ @@ -10338,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*/ @@ -10354,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; @@ -10368,35 +10591,37 @@ 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 */ }/* End loop cptcovt */ /* Checking for missing or useless values in comparison of current model needs */ /* Feeds resultmodel[nres][k1]=k2 for single covariate (k1) in the model equation */ - for(k2=1; k2 <=j;k2++){ /* Loop on resultline variables: result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */ + for(k2=1; k2 <=j;k2++){ /* j or cptcovs is the number of single covariates used either in the model line as well as in the result line (dummy or quantitative) + * Loop on resultline variables: result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */ match=0; for(k1=1; k1<= cptcovt ;k1++){ /* loop on model: model line V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ if(Typevar[k1]==0){ /* Single only */ if(Tvar[k1]==Tvarsel[k2]) { /* Tvar[2]=4 == Tvarsel[1]=4 */ resultmodel[nres][k1]=k2; /* k1th position in the model equation corresponds to k2th position in the result line. resultmodel[2]=1 resultmodel[1]=2 resultmodel[3]=3 resultmodel[6]=4 resultmodel[9]=5 */ + modelresult[nres][k2]=k1; /* k1th position in the model equation corresponds to k2th position in the result line. modelresult[1]=2 modelresult[2]=1 modelresult[3]=3 remodelresult[4]=6 modelresult[5]=9 */ ++match; } } } 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; } } - + /* cptcovres=j /\* Number of variables in the resultline is equal to cptcovs and thus useless *\/ */ /* We need to deduce which combination number is chosen and save quantitative values */ /* model line V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* nres=1st result line: V4=1 V5=25.1 V3=0 V2=8 V1=1 */ @@ -10415,12 +10640,13 @@ int decoderesult( char resultline[], int /* V(Tvqresult)=Tqresult V5=25.1 V2=8 Tqresult[nres=1][1]=25.1 */ /* V5*age V5 known which value for nres? */ /* Tqinvresult[2]=8 Tqinvresult[1]=25.1 */ - for(k1=1, k=0, k4=0, k4q=0; k1 <=cptcovt;k1++){ /* loop k1 on position in the model line (excluding product) */ + for(k1=1, k=0, k4=0, k4q=0; k1 <=cptcovt;k1++){ /* cptcovt number of covariates (excluding 1 and age or age*age) in the MODEL equation. + * loop on position k1 in the MODEL LINE */ /* k counting number of combination of single dummies in the equation model */ /* k4 counting single dummies in the equation model */ /* k4q counting single quantitatives in the equation model */ - if( Dummy[k1]==0 && Typevar[k1]==0 ){ /* Dummy and Single */ - /* k4+1= position in the resultline V(Tvarsel)=Tvalsel=Tresult[nres][pos](value); V(Tvresult[nres][pos] (variable): V(variable)=value) */ + if( Dummy[k1]==0 && Typevar[k1]==0 ){ /* Dummy and Single, k1 is sorting according to MODEL, but k3 to resultline */ + /* k4+1= (not always if quant in model) position in the resultline V(Tvarsel)=Tvalsel=Tresult[nres][pos](value); V(Tvresult[nres][pos] (variable): V(variable)=value) */ /* modelresult[k3]=k1: k3th position in the result line corresponds to the k1 position in the model line (doesn't work with products)*/ /* Value in the (current nres) resultline of the variable at the k1th position in the model equation resultmodel[nres][k1]= k3 */ /* resultmodel[nres][k1]=k3: k1th position in the model correspond to the k3 position in the resultline */ @@ -10428,19 +10654,21 @@ int decoderesult( char resultline[], int /* Tvarsel[k3]: Name of the variable at the k3th position in the result line. */ /* Tvalsel[k3]: Value of the variable at the k3th position in the result line. */ /* Tresult[nres][result_position]= value of the dummy variable at the result_position in the nres resultline */ - /* Tvresult[nres][result_position]= id of the dummy variable at the result_position in the nres resultline */ + /* Tvresult[nres][result_position]= name of the dummy variable at the result_position in the nres resultline */ /* Tinvresult[nres][Name of a dummy variable]= value of the variable in the result line */ /* TinvDoQresult[nres][Name of a Dummy or Q variable]= value of the variable in the result line */ k3= resultmodel[nres][k1]; /* From position k1 in model get position k3 in result line */ /* nres=1 k1=2 resultmodel[2(V4)] = 1=k3 ; k1=3 resultmodel[3(V3)] = 2=k3*/ k2=(int)Tvarsel[k3]; /* from position k3 in resultline get name k2: nres=1 k1=2=>k3=1 Tvarsel[resultmodel[2]]= Tvarsel[1] = 4=k2 (V4); k1=3=>k3=2 Tvarsel[2]=3 (V3)*/ k+=Tvalsel[k3]*pow(2,k4); /* nres=1 k1=2 Tvalsel[1]=1 (V4=1); k1=3 k3=2 Tvalsel[2]=0 (V3=0) */ - TinvDoQresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* Stores the value into the name of the variable. */ + TinvDoQresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* TinvDoQresult[nres][Name]=Value; stores the value into the name of the variable. */ /* Tinvresult[nres][4]=1 */ - Tresult[nres][k4+1]=Tvalsel[k3];/* Tresult[nres=2][1]=1(V4=1) Tresult[nres=2][2]=0(V3=0) */ - Tvresult[nres][k4+1]=(int)Tvarsel[k3];/* Tvresult[nres][1]=4 Tvresult[nres][3]=1 */ + /* Tresult[nres][k4+1]=Tvalsel[k3];/\* Tresult[nres=2][1]=1(V4=1) Tresult[nres=2][2]=0(V3=0) *\/ */ + Tresult[nres][k3]=Tvalsel[k3];/* Tresult[nres=2][1]=1(V4=1) Tresult[nres=2][2]=0(V3=0) */ + /* Tvresult[nres][k4+1]=(int)Tvarsel[k3];/\* Tvresult[nres][1]=4 Tvresult[nres][3]=1 *\/ */ + Tvresult[nres][k3]=(int)Tvarsel[k3];/* Tvresult[nres][1]=4 Tvresult[nres][3]=1 */ Tinvresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* Tinvresult[nres][4]=1 */ - precov[nres][k1]=Tvalsel[k3]; + precov[nres][k1]=Tvalsel[k3]; /* Value from resultline of the variable at the k1 position in the model */ printf("Decoderesult Dummy k=%d, k1=%d precov[nres=%d][k1=%d]=%.f V(k2=V%d)= Tvalsel[%d]=%d, 2**(%d)\n",k, k1, nres, k1,precov[nres][k1], k2, k3, (int)Tvalsel[k3], k4); k4++;; }else if( Dummy[k1]==1 && Typevar[k1]==0 ){ /* Quantitative and single */ @@ -10450,8 +10678,12 @@ int decoderesult( char resultline[], int k3q= resultmodel[nres][k1]; /* resultmodel[1(V5)] = 5 =k3q */ k2q=(int)Tvarsel[k3q]; /* Name of variable at k3q th position in the resultline */ /* Tvarsel[resultmodel[1]]= Tvarsel[1] = 4=k2 */ - Tqresult[nres][k4q+1]=Tvalsel[k3q]; /* Tqresult[nres][1]=25.1 */ - Tvqresult[nres][k4q+1]=(int)Tvarsel[k3q]; /* Tvqresult[nres][1]=5 */ + /* Tqresult[nres][k4q+1]=Tvalsel[k3q]; /\* Tqresult[nres][1]=25.1 *\/ */ + /* Tvresult[nres][k4q+1]=(int)Tvarsel[k3q];/\* Tvresult[nres][1]=4 Tvresult[nres][3]=1 *\/ */ + /* Tvqresult[nres][k4q+1]=(int)Tvarsel[k3q]; /\* Tvqresult[nres][1]=5 *\/ */ + Tqresult[nres][k3q]=Tvalsel[k3q]; /* Tqresult[nres][1]=25.1 */ + Tvresult[nres][k3q]=(int)Tvarsel[k3q];/* Tvresult[nres][1]=4 Tvresult[nres][3]=1 */ + Tvqresult[nres][k3q]=(int)Tvarsel[k3q]; /* Tvqresult[nres][1]=5 */ Tqinvresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* Tqinvresult[nres][5]=25.1 */ TinvDoQresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* Tqinvresult[nres][5]=25.1 */ precov[nres][k1]=Tvalsel[k3q]; @@ -10463,15 +10695,15 @@ int decoderesult( char resultline[], int k3= resultmodel[nres][k1]; /* nres=1 k1=2 resultmodel[2(V4)] = 1=k3 ; k1=3 resultmodel[3(V3)] = 2=k3*/ k2=(int)Tvarsel[k3]; /* nres=1 k1=2=>k3=1 Tvarsel[resultmodel[2]]= Tvarsel[1] = 4=k2 (V4); k1=3=>k3=2 Tvarsel[2]=3 (V3)*/ - TinvDoQresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* Tinvresult[nres][4]=1 */ + TinvDoQresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* TinvDoQresult[nres][4]=1 */ precov[nres][k1]=Tvalsel[k3]; printf("Decoderesult Dummy with age k=%d, k1=%d precov[nres=%d][k1=%d]=%.f Tvar[%d]=V%d k2=Tvarsel[%d]=%d Tvalsel[%d]=%d\n",k, k1, nres, k1,precov[nres][k1], k1, Tvar[k1], k3,(int)Tvarsel[k3], k3, (int)Tvalsel[k3]); }else if( Dummy[k1]==3 ){ /* For quant with age product */ k3q= resultmodel[nres][k1]; /* resultmodel[1(V5)] = 25.1=k3q */ k2q=(int)Tvarsel[k3q]; /* Tvarsel[resultmodel[1]]= Tvarsel[1] = 4=k2 */ - TinvDoQresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* Tqinvresult[nres][5]=25.1 */ + TinvDoQresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* TinvDoQresult[nres][5]=25.1 */ precov[nres][k1]=Tvalsel[k3q]; - printf("Decoderesult Quantitative with age nres=%d, k1=%d, precov[nres=%d][k1=%d]=%.f Tvar[%d]=V%d V(k2q=%d)= Tvarsel[%d]=%d, Tvalsel[%d]=%f\n",nres, k1, nres, k1,precov[nres][k1], k1, Tvar[k1], k2q, k3q, Tvarsel[k3q], k3q, Tvalsel[k3q]); + printf("Decoderesult Quantitative with age nres=%d, k1=%d, precov[nres=%d][k1=%d]=%f Tvar[%d]=V%d V(k2q=%d)= Tvarsel[%d]=%d, Tvalsel[%d]=%f\n",nres, k1, nres, k1,precov[nres][k1], k1, Tvar[k1], k2q, k3q, Tvarsel[k3q], k3q, Tvalsel[k3q]); }else if(Typevar[k1]==2 ){ /* For product quant or dummy (not with age) */ precov[nres][k1]=TinvDoQresult[nres][Tvardk[k1][1]] * TinvDoQresult[nres][Tvardk[k1][2]]; printf("Decoderesult Quantitative or Dummy (not with age) nres=%d k1=%d precov[nres=%d][k1=%d]=%.f V%d(=%.f) * V%d(=%.f) \n",nres, k1, nres, k1,precov[nres][k1], Tvardk[k1][1], TinvDoQresult[nres][Tvardk[k1][1]], Tvardk[k1][2], TinvDoQresult[nres][Tvardk[k1][2]]); @@ -10481,7 +10713,7 @@ int decoderesult( char resultline[], int } } - TKresult[nres]=++k; /* Combination for the nresult and the model */ + TKresult[nres]=++k; /* Number of combinations of dummies for the nresult and the model =Tvalsel[k3]*pow(2,k4) + 1*/ return (0); } @@ -10517,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; @@ -10641,9 +10873,17 @@ int decodemodel( char model[], int lasto Tvar[k]=ncovcol+nqv+ntv+nqtv+k1; /* For model-covariate k tells which data-covariate to use but because this model-covariate is a construction we invent a new column which is after existing variables ncovcol+nqv+ntv+nqtv + k1 - If already ncovcol=4 and model=V2 + V1 +V1*V4 +age*V3 +V3*V2 + If already ncovcol=4 and model= V2 + V1 + V1*V4 + age*V3 + V3*V2 thus after V4 we invent V5 and V6 because age*V3 will be computed in 4 Tvar[3=V1*V4]=4+1=5 Tvar[5=V3*V2]=4 + 2= 6, Tvar[4=age*V3]=4 etc */ + /* Please remark that the new variables are model dependent */ + /* If we have 4 variable but the model uses only 3, like in + * model= V1 + age*V1 + V2 + V3 + age*V2 + age*V3 + V1*V2 + V1*V3 + * k= 1 2 3 4 5 6 7 8 + * Tvar[k]=1 1 2 3 2 3 (5 6) (and not 4 5 because of V4 missing) + * Tage[kk] [1]= 2 [2]=5 [3]=6 kk=1 to cptcovage=3 + * Tvar[Tage[kk]][1]=2 [2]=2 [3]=3 + */ Typevar[k]=2; /* 2 for double fixed dummy covariates */ cutl(strc,strb,strd,'V'); /* strd was Vm, strc is m */ Tprod[k1]=k; /* Tprod[1]=3(=V1*V4) for V2+V1+V1*V4+age*V3+V3*V2 */ @@ -10746,8 +10986,8 @@ Dummy[k] 0=dummy (0 1), 1 quantitative ( modell[k].maintype= FTYPE; modell[k].subtype= FQ; nsq++; - TvarsQ[nsq]=Tvar[k]; - TvarsQind[nsq]=k; + TvarsQ[nsq]=Tvar[k]; /* Gives the variable name (extended to products) of first nsq simple quantitative covariates (fixed or time vary see below */ + TvarsQind[nsq]=k; /* Gives the position in the model equation of the first nsq simple quantitative covariates (fixed or time vary) */ ncovf++; TvarF[ncovf]=Tvar[k]; TvarFind[ncovf]=k; @@ -10778,8 +11018,8 @@ Dummy[k] 0=dummy (0 1), 1 quantitative ( modell[k].subtype= VQ; ncovv++; /* Only simple time varying variables */ nsq++; - TvarsQ[nsq]=Tvar[k]; /* k=1 Tvar=5 nsq=1 TvarsQ[1]=5 */ - TvarsQind[nsq]=k; /* For single quantitative covariate gives the model position of each single quantitative covariate */ + TvarsQ[nsq]=Tvar[k]; /* k=1 Tvar=5 nsq=1 TvarsQ[1]=5 */ /* Gives the variable name (extended to products) of first nsq simple quantitative covariates (fixed or time vary here) */ + TvarsQind[nsq]=k; /* For single quantitative covariate gives the model position of each single quantitative covariate *//* Gives the position in the model equation of the first nsq simple quantitative covariates (fixed or time vary) */ TvarV[ncovv]=Tvar[k]; TvarVind[ncovv]=k; /* TvarVind[1]=1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Any time varying singele */ TvarVQ[nqtveff]=Tvar[k]; /* TvarVQ[1]=V5 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple time varying quantitative variable */ @@ -10969,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); } } @@ -11334,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++){ */ @@ -11348,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"); @@ -11370,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]; @@ -11390,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; } @@ -11433,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"); @@ -11461,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"); @@ -11485,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]; @@ -11496,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); */ @@ -11506,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; @@ -11516,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; @@ -11618,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,"\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," %.5f", p3mat[i][j][h]);/* Bug valgrind */ 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 */ @@ -11722,7 +11986,7 @@ int main(int argc, char *argv[]) char pathr[MAXLINE], pathimach[MAXLINE]; char *tok, *val; /* pathtot */ - int firstobs=1, lastobs=10; /* nobs = lastobs-firstobs declared globally ;*/ + /* int firstobs=1, lastobs=10; /\* nobs = lastobs-firstobs declared globally ;*\/ */ int c, h , cpt, c2; int jl=0; int i1, j1, jk, stepsize=0; @@ -12061,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); @@ -12345,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 */ @@ -12365,6 +12629,7 @@ Please run with mle=-1 to get a correct Tvar=ivector(1,NCOVMAX); /* Was 15 changed to NCOVMAX. */ TvarsDind=ivector(1,NCOVMAX); /* */ TnsdVar=ivector(1,NCOVMAX); /* */ + /* for(i=1; i<=NCOVMAX;i++) TnsdVar[i]=3; */ TvarsD=ivector(1,NCOVMAX); /* */ TvarsQind=ivector(1,NCOVMAX); /* */ TvarsQ=ivector(1,NCOVMAX); /* */ @@ -12487,12 +12752,12 @@ Please run with mle=-1 to get a correct Ndum =ivector(-1,NCOVMAX); cptcoveff=0; if (ncovmodel-nagesqr > 2 ){ /* That is if covariate other than cst, age and age*age */ - tricode(&cptcoveff,Tvar,nbcode,imx, Ndum); /**< Fills nbcode[Tvar[j]][l]; */ + tricode(&cptcoveff,Tvar,nbcode,imx, Ndum); /**< Fills nbcode[Tvar[j]][l]; as well as calculate cptcoveff or number of total effective dummy covariates*/ } ncovcombmax=pow(2,cptcoveff); - invalidvarcomb=ivector(1, ncovcombmax); - for(i=1;iDatafile=%s Firstpass=%d La optionfilehtmcov,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model); } - fprintf(fichtm,"\n\n\nIMaCh %s\n IMaCh for Interpolated Markov Chain
\nSponsored by Copyright (C) 2002-2015 INED-EUROREVES-Institut de longévité-2013-2022-Japan Society for the Promotion of Sciences 日本学術振興会 (Grant-in-Aid for Scientific Research 25293121) - Intel Software 2015-2018
\ -
\n\ + fprintf(fichtm,"\n\n\ +IMaCh %s\n\ + IMaCh for Interpolated Markov Chain
\n\ +Sponsored by Copyright (C) 2002-2015 INED\ +-EUROREVES-Institut de longévité-2013-2022-Japan Society for the Promotion of Sciences 日本学術振興会 \ +(Grant-in-Aid for Scientific Research 25293121) - \ +Intel Software 2015-2018
\n", optionfilehtm); + + fprintf(fichtm,"
\n\ IMaCh-%s
%s \
\n\ -Title=%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\ @@ -12639,7 +12911,7 @@ Title=%s
    Datafile=%s Firstpass=%d La - Log file of the run: %s
    \n\ - Gnuplot file name: %s
    \n\ - Date and time at start: %s
\n",\ - optionfilehtm,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model,\ + version,fullversion,optionfilehtm,optionfilehtm,title,datafile,datafile,firstpass,lastpass,stepm, weightopt, model, \ optionfilefiname,optionfilext,optionfilefiname,optionfilext,\ fileres,fileres,\ filelog,filelog,optionfilegnuplot,optionfilegnuplot,strstart); @@ -12953,6 +13225,7 @@ Please run with mle=-1 to get a correct globpr=1; /* again, to print the individual contributions using computed gpimx and gsw */ likelione(ficres, p, npar, nlstate, &globpr, &ipmx, &sw, &fretone, funcone); /* Prints the contributions to the likelihood */ printf("Second Likeli=%12.6f ipmx=%ld sw=%12.6f",fretone,ipmx,sw); + /* exit(0); */ for (k=1; k<=npar;k++) printf(" %d %8.5f",k,p[k]); printf("\n"); @@ -13631,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 ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficreseij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][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"); @@ -13650,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) { @@ -13689,36 +13962,88 @@ Please run with mle=-1 to get a correct i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */ if (cptcovn < 1){i1=1;} - for(nres=1; nres <= nresult; nres++) /* For each resultline */ - for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */ - if(i1 != 1 && TKresult[nres]!= k) + for(nres=1; nres <= nresult; nres++) /* For each resultline, find the combination and output results according to the values of dummies and then quanti. */ + for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying. For each nres and each value at position k + * we know Tresult[nres][result_position]= value of the dummy variable at the result_position in the nres resultline + * Tvqresult[nres][result_position]= id of the variable at the result_position in the nres resultline + * and Tqresult[nres][result_position]= value of the variable at the result_position in the nres resultline */ + /* */ + if(i1 != 1 && TKresult[nres]!= k) /* TKresult[nres] is the combination of this nres resultline. All the i1 combinations are not output */ continue; printf("\n# model %s \n#****** Result for:", model); fprintf(ficrest,"\n# model %s \n#****** Result for:", model); fprintf(ficlog,"\n# model %s \n#****** Result for:", model); - for(j=1;j<=cptcoveff;j++){ - printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); - fprintf(ficrest,"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(ficrest," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - } + /* It might not be a good idea to mix dummies and quantitative */ + /* for(j=1;j<=cptcoveff;j++){ /\* j=resultpos. Could be a loop on cptcovs: number of single dummy covariate in the result line as well as in the model *\/ */ + for(j=1;j<=cptcovs;j++){ /* j=resultpos. Could be a loop on cptcovs: number of single covariate (dummy or quantitative) in the result line as well as in the model */ + /* printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); /\* Output by variables in the resultline *\/ */ + /* Tvaraff[j] is the name of the dummy variable in position j in the equation model: + * Tvaraff[1]@9={4, 3, 0, 0, 0, 0, 0, 0, 0}, in model=V5+V4+V3+V4*V3+V5*age + * (V5 is quanti) V4 and V3 are dummies + * TnsdVar[4] is the position 1 and TnsdVar[3]=2 in codtabm(k,l)(V4 V3)=V4 V3 + * l=1 l=2 + * k=1 1 1 0 0 + * k=2 2 1 1 0 + * k=3 [1] [2] 0 1 + * k=4 2 2 1 1 + * If nres=1 result: V3=1 V4=0 then k=3 and outputs + * If nres=2 result: V4=1 V3=0 then k=2 and outputs + * nres=1 =>k=3 j=1 V4= nbcode[4][codtabm(3,1)=1)=0; j=2 V3= nbcode[3][codtabm(3,2)=2]=1 + * nres=2 =>k=2 j=1 V4= nbcode[4][codtabm(2,1)=2)=1; j=2 V3= nbcode[3][codtabm(2,2)=1]=0 + */ + /* Tvresult[nres][j] Name of the variable at position j in this resultline */ + /* Tresult[nres][j] Value of this variable at position j could be a float if quantitative */ +/* We give up with the combinations!! */ + printf("\n j=%d In computing T_ Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=%d cptcovs=%d, cptcoveff=%d Fixed[modelresult[nres][j]]=%d\n", j, nres, j, Dummy[modelresult[nres][j]],nres,j,modelresult[nres][j],cptcovs, cptcoveff,Fixed[modelresult[nres][j]]); /* end if dummy or quanti */ + + if(Dummy[modelresult[nres][j]]==0){/* Dummy variable of the variable in position modelresult in the model corresponding to j in resultline */ + printf("V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + fprintf(ficlog,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + fprintf(ficrest,"V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][j]); /* Output of each value for the combination TKresult[nres], ordere by the covariate values in the resultline */ + if(Fixed[modelresult[nres][j]]==0){ /* Fixed */ + printf("fixed ");fprintf(ficlog,"fixed ");fprintf(ficrest,"fixed "); + }else{ + printf("varyi ");fprintf(ficlog,"varyi ");fprintf(ficrest,"varyi "); + } + /* fprintf(ficrest,"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]])]); */ + }else if(Dummy[modelresult[nres][j]]==1){ /* Quanti variable */ + /* For each selected (single) quantitative value */ + 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{ + printf("varyi ");fprintf(ficlog,"varyi ");fprintf(ficrest,"varyi "); + } + }else{ + printf("Error in computing T_ Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=%d cptcovs=%d, cptcoveff=%d \n", nres, j, Dummy[modelresult[nres][j]],nres,j,modelresult[nres][j],cptcovs, cptcoveff); /* end if dummy or quanti */ + fprintf(ficlog,"Error in computing T_ Dummy[modelresult[%d][%d]]=%d, modelresult[%d][%d]=%d cptcovs=%d, cptcoveff=%d \n", nres, j, Dummy[modelresult[nres][j]],nres,j,modelresult[nres][j],cptcovs, cptcoveff); /* end if dummy or quanti */ + exit(1); + } + } /* End loop for each variable in the resultline */ + /* for (j=1; j<= nsq; j++){ /\* For each selected (single) quantitative value *\/ */ + /* printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); /\* Wrong j is not in the equation model *\/ */ + /* fprintf(ficrest," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); */ + /* } */ fprintf(ficrest,"******\n"); fprintf(ficlog,"******\n"); printf("******\n"); 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 ",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 */ - fprintf(ficresstdeij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); - fprintf(ficrescveij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][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=%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"); @@ -13726,10 +14051,11 @@ Please run with mle=-1 to get a correct fprintf(ficresvij,"\n#****** "); /* pstamp(ficresvij); */ for(j=1;j<=cptcoveff;j++) - fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[TnsdVar[Tvaraff[j]]])]); + fprintf(ficresvij,"V%d=%d ",Tvresult[nres][j],Tresult[nres][j]); + /* fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[TnsdVar[Tvaraff[j]]])]); */ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */ /* fprintf(ficresvij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]); /\* To solve *\/ */ - fprintf(ficresvij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); /* Solved */ + fprintf(ficresvij," V%d=%lg ",Tvqresult[nres][j],Tqresult[nres][resultmodel[nres][j]]); /* Solved */ } fprintf(ficresvij,"******\n"); @@ -13760,7 +14086,7 @@ Please run with mle=-1 to get a correct fprintf(ficrest,"the age specific prevalence observed (cross-sectionally) in the population i.e cross-sectionally\n in each health state (popbased=1) (mobilav=%d)\n",mobilav); else fprintf(ficrest,"the age specific forward period (stable) prevalences in each health state \n"); - fprintf(ficrest,"# Age popbased mobilav e.. (std) "); + fprintf(ficrest,"# Age popbased mobilav e.. (std) "); /* Adding covariate values? */ for (i=1;i<=nlstate;i++) fprintf(ficrest,"e.%d (std) ",i); fprintf(ficrest,"\n"); /* printf("Which p?\n"); for(i=1;i<=npar;i++)printf("p[i=%d]=%lf,",i,p[i]);printf("\n"); */ @@ -13807,12 +14133,12 @@ Please run with mle=-1 to get a correct printf("done selection\n");fflush(stdout); fprintf(ficlog,"done selection\n");fflush(ficlog); - } /* End k selection */ + } /* End k selection or end covariate selection for nres */ printf("done State-specific expectancies\n");fflush(stdout); fprintf(ficlog,"done State-specific expectancies\n");fflush(ficlog); - /* variance-covariance of forward period prevalence*/ + /* variance-covariance of forward period prevalence */ varprlim(fileresu, nresult, mobaverage, mobilavproj, bage, fage, prlim, &ncvyear, ftolpl, p, matcov, delti, stepm, cptcoveff); @@ -13886,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);