--- imach/src/imach.c 2016/08/23 16:51:20 1.234
+++ imach/src/imach.c 2016/08/26 09:20:19 1.237
@@ -1,6 +1,15 @@
-/* $Id: imach.c,v 1.234 2016/08/23 16:51:20 brouard Exp $
+/* $Id: imach.c,v 1.237 2016/08/26 09:20:19 brouard Exp $
$State: Exp $
$Log: imach.c,v $
+ Revision 1.237 2016/08/26 09:20:19 brouard
+ Summary: to valgrind
+
+ Revision 1.236 2016/08/25 10:50:18 brouard
+ *** empty log message ***
+
+ Revision 1.235 2016/08/25 06:59:23 brouard
+ *** empty log message ***
+
Revision 1.234 2016/08/23 16:51:20 brouard
*** empty log message ***
@@ -902,12 +911,12 @@ typedef struct {
#define ODIRSEPARATOR '\\'
#endif
-/* $Id: imach.c,v 1.234 2016/08/23 16:51:20 brouard Exp $ */
+/* $Id: imach.c,v 1.237 2016/08/26 09:20:19 brouard Exp $ */
/* $State: Exp $ */
#include "version.h"
char version[]=__IMACH_VERSION__;
char copyright[]="February 2016,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2018";
-char fullversion[]="$Revision: 1.234 $ $Date: 2016/08/23 16:51:20 $";
+char fullversion[]="$Revision: 1.237 $ $Date: 2016/08/26 09:20:19 $";
char strstart[80];
char optionfilext[10], optionfilefiname[FILENAMELENGTH];
int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */
@@ -1105,6 +1114,16 @@ int *TvarsDind;
int *TvarsQ;
int *TvarsQind;
+#define MAXRESULTLINES 10
+int nresult=0;
+int TKresult[MAXRESULTLINES];
+int Tresult[MAXRESULTLINES][NCOVMAX];/* For dummy variable , value (output) */
+int Tinvresult[MAXRESULTLINES][NCOVMAX];/* For dummy variable , value (output) */
+int Tvresult[MAXRESULTLINES][NCOVMAX]; /* For dummy variable , variable # (output) */
+double Tqresult[MAXRESULTLINES][NCOVMAX]; /* For quantitative variable , value (output) */
+double Tqinvresult[MAXRESULTLINES][NCOVMAX]; /* For quantitative variable , value (output) */
+int Tvqresult[MAXRESULTLINES][NCOVMAX]; /* For quantitative variable , variable # (output) */
+
/* int *TDvar; /\**< TDvar[1]=4, TDvarF[2]=3, TDvar[3]=6 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 *\/ */
int *TvarF; /**< TvarF[1]=Tvar[6]=2, TvarF[2]=Tvar[7]=7, TvarF[3]=Tvar[9]=1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
int *TvarFind; /**< TvarFind[1]=6, TvarFind[2]=7, Tvarind[3]=9 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
@@ -2333,9 +2352,11 @@ void powell(double p[], double **xi, int
/**** Prevalence limit (stable or period prevalence) ****************/
-double **prevalim(double **prlim, int nlstate, double x[], double age, double **oldm, double **savm, double ftolpl, int *ncvyear, int ij)
+ 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 combiation ij by left multiplying the unit
+ /* Computes the prevalence limit in each live state at age x and for covariate combination ij
+ (and selected quantitative values in nres)
+ by left multiplying the unit
matrix by transitions matrix until convergence is reached with precision ftolpl */
/* Wx= Wx-1 Px-1= Wx-2 Px-2 Px-1 = Wx-n Px-n ... Px-2 Px-1 I */
/* Wx is row vector: population in state 1, population in state 2, population dead */
@@ -2387,27 +2408,36 @@ double **prevalim(double **prlim, int nl
for (k=1; k<=nsd;k++) { /* For single dummy covariates only */
/* Here comes the value of the covariate 'ij' after renumbering k with single dummy covariates */
cov[2+nagesqr+TvarsDind[k]]=nbcode[TvarsD[k]][codtabm(ij,k)];
- printf("prevalim ij=%d k=%d TvarsD[%d]=%d TvarsDind[%d]=%d nbcode=%d cov=%lf codtabm(%d,Tvar[%d])=%d \n",ij,k, k, TvarsD[k],k,TvarsDind[k],nbcode[TvarsD[k]][codtabm(ij,k)],cov[2+nagesqr+TvarsDind[k]], ij, k, codtabm(ij,k));
+ /* printf("prevalim Dummy combi=%d k=%d TvarsD[%d]=V%d TvarsDind[%d]=%d nbcode=%d cov=%lf codtabm(%d,Tvar[%d])=%d \n",ij,k, k, TvarsD[k],k,TvarsDind[k],nbcode[TvarsD[k]][codtabm(ij,k)],cov[2+nagesqr+TvarsDind[k]], ij, k, codtabm(ij,k)); */
}
for (k=1; k<=nsq;k++) { /* For single varying covariates only */
/* Here comes the value of quantitative after renumbering k with single quantitative covariates */
- /* cov[2+nagesqr+TvarsQind[k]]=qselvar[k]; */
- printf("prevalim ij=%d k=%d TvarsQind[%d]=%d \n",ij,k,k,TvarsQind[k]);
+ cov[2+nagesqr+TvarsQind[k]]=Tqresult[nres][k];
+ /* printf("prevalim Quantitative k=%d TvarsQind[%d]=%d, TvarsQ[%d]=V%d,Tqresult[%d][%d]=%f\n",k,k,TvarsQind[k],k,TvarsQ[k],nres,k,Tqresult[nres][k]); */
}
- /*wrong? for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */
- /* for (k=1; k<=cptcovage;k++) cov[2+nagesqr+Tage[k]]=nbcode[Tvar[k]][codtabm(ij,Tvar[k])]*cov[2]; */
- for (k=1; k<=cptcovage;k++){
+ for (k=1; k<=cptcovage;k++){ /* For product with age */
if(Dummy[Tvar[Tage[k]]]){
cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2];
} else{
- ;
- /* cov[2+nagesqr+Tage[k]]=qselvar[k]; */
+ cov[2+nagesqr+Tage[k]]=Tqresult[nres][k];
}
- printf("prevalim Age ij=%d k=%d Tage[%d]=%d \n",ij,k,k,Tage[k]);
+ /* printf("prevalim Age combi=%d k=%d Tage[%d]=V%d Tqresult[%d][%d]=%f\n",ij,k,k,Tage[k],nres,k,Tqresult[nres][k]); */
}
- for (k=1; k<=cptcovprod;k++){ /* */
- printf("prevalim Prod ij=%d k=%d Tprod[%d]=%d Tvard[%d][1]=%d, Tvard[%d][2]=%d\n",ij,k,k,Tprod[k], k,Tvard[k][1], k,Tvard[k][2]);
- cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)] * nbcode[Tvard[k][2]][codtabm(ij,k)];
+ for (k=1; k<=cptcovprod;k++){ /* For product without age */
+ /* printf("prevalim Prod ij=%d k=%d Tprod[%d]=%d Tvard[%d][1]=V%d, Tvard[%d][2]=V%d\n",ij,k,k,Tprod[k], k,Tvard[k][1], k,Tvard[k][2]); */
+ if(Dummy[Tvard[k][1]==0]){
+ if(Dummy[Tvard[k][2]==0]){
+ cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)] * nbcode[Tvard[k][2]][codtabm(ij,k)];
+ }else{
+ cov[2+nagesqr+Tprod[k]]=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(ij,k)] * Tqinvresult[nres][Tvard[k][1]];
+ }else{
+ cov[2+nagesqr+Tprod[k]]=Tqinvresult[nres][Tvard[k][1]]* Tqinvresult[nres][Tvard[k][2]];
+ }
+ }
}
/*printf("ij=%d cptcovprod=%d tvar=%d ", ij, cptcovprod, Tvar[1]);*/
/*printf("ij=%d cov[3]=%lf cov[4]=%lf \n",ij, cov[3],cov[4]);*/
@@ -2530,8 +2560,6 @@ Earliest age to start was %d-%d=%d, ncvl
cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,k)];
/* printf("prevalim ij=%d k=%d Tvar[%d]=%d nbcode=%d cov=%lf codtabm(%d,Tvar[%d])=%d \n",ij,k, k, Tvar[k],nbcode[Tvar[k]][codtabm(ij,Tvar[k])],cov[2+k], ij, k, codtabm(ij,Tvar[k])]); */
}
- /*wrong? for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */
- /* for (k=1; k<=cptcovage;k++) cov[2+nagesqr+Tage[k]]=nbcode[Tvar[k]][codtabm(ij,Tvar[k])]*cov[2]; */
for (k=1; k<=cptcovage;k++) cov[2+nagesqr+Tage[k]]=nbcode[Tvar[k]][codtabm(ij,k)]*cov[2];
for (k=1; k<=cptcovprod;k++) /* Useless */
/* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,Tvard[k][1])] * nbcode[Tvard[k][2]][codtabm(ij,Tvard[k][2])]; */
@@ -2856,7 +2884,7 @@ double **matprod2(double **out, double *
/************* Higher Matrix Product ***************/
-double ***hpxij(double ***po, int nhstepm, double age, int hstepm, double *x, int nlstate, int stepm, double **oldm, double **savm, int ij )
+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 combination of covariate values corresponding to ij over
'nhstepm*hstepm*stepm' months (i.e. until
@@ -2892,16 +2920,34 @@ double ***hpxij(double ***po, int nhstep
cov[2]=agexact;
if(nagesqr==1)
cov[3]= agexact*agexact;
- for (k=1; k<=cptcovn;k++)
- cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,k)];
- /* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,Tvar[k])]; */
- for (k=1; k<=cptcovage;k++) /* Should start at cptcovn+1 */
- /* cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */
- cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2];
- /* cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,Tvar[Tage[k]])]*cov[2]; */
- for (k=1; k<=cptcovprod;k++) /* Useless because included in cptcovn */
- cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)];
- /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,Tvard[k][1])]*nbcode[Tvard[k][2]][codtabm(ij,Tvard[k][2])]; */
+ for (k=1; k<=nsd;k++) { /* For single dummy covariates only */
+ /* Here comes the value of the covariate 'ij' after renumbering k with single dummy covariates */
+ cov[2+nagesqr+TvarsDind[k]]=nbcode[TvarsD[k]][codtabm(ij,k)];
+ /* printf("hpxij Dummy combi=%d k=%d TvarsD[%d]=V%d TvarsDind[%d]=%d nbcode=%d cov=%lf codtabm(%d,Tvar[%d])=%d \n",ij,k, k, TvarsD[k],k,TvarsDind[k],nbcode[TvarsD[k]][codtabm(ij,k)],cov[2+nagesqr+TvarsDind[k]], ij, k, codtabm(ij,k)); */
+ }
+ for (k=1; k<=nsq;k++) { /* For single varying covariates only */
+ /* Here comes the value of quantitative after renumbering k with single quantitative covariates */
+ cov[2+nagesqr+TvarsQind[k]]=Tqresult[nres][k];
+ /* printf("hPxij Quantitative k=%d TvarsQind[%d]=%d, TvarsQ[%d]=V%d,Tqresult[%d][%d]=%f\n",k,k,TvarsQind[k],k,TvarsQ[k],nres,k,Tqresult[nres][k]); */
+ }
+ for (k=1; k<=cptcovage;k++){
+ if(Dummy[Tvar[Tage[k]]]){
+ cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2];
+ } else{
+ cov[2+nagesqr+Tage[k]]=Tqresult[nres][k];
+ }
+ /* printf("hPxij Age combi=%d k=%d Tage[%d]=V%d Tqresult[%d][%d]=%f\n",ij,k,k,Tage[k],nres,k,Tqresult[nres][k]); */
+ }
+ for (k=1; k<=cptcovprod;k++){ /* */
+ /* printf("hPxij Prod ij=%d k=%d Tprod[%d]=%d Tvard[%d][1]=V%d, Tvard[%d][2]=V%d\n",ij,k,k,Tprod[k], k,Tvard[k][1], k,Tvard[k][2]); */
+ cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)] * nbcode[Tvard[k][2]][codtabm(ij,k)];
+ }
+ /* for (k=1; k<=cptcovn;k++) */
+ /* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,k)]; */
+ /* for (k=1; k<=cptcovage;k++) /\* Should start at cptcovn+1 *\/ */
+ /* cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; */
+ /* for (k=1; k<=cptcovprod;k++) /\* Useless because included in cptcovn *\/ */
+ /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)]; */
/*printf("hxi cptcov=%d cptcode=%d\n",cptcov,cptcode);*/
@@ -4058,7 +4104,7 @@ void freqsummary(char fileres[], int ia
Title=%s
Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s
\n",\
fileresphtm,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model);
}
- fprintf(ficresphtm,"Current page is file %s
\n\n
Frequencies and prevalence by age at begin of transition
\n",fileresphtm, fileresphtm);
+ fprintf(ficresphtm,"Current page is file %s
\n\nFrequencies and prevalence by age at begin of transition and dummy covariate value at beginning of transition
\n",fileresphtm, fileresphtm);
strcpy(fileresphtmfr,subdirfext(optionfilefiname,"PHTMFR_",".htm"));
if((ficresphtmfr=fopen(fileresphtmfr,"w"))==NULL) {
@@ -4862,7 +4908,7 @@ void concatwav(int wav[], int **dh, int
/*********** Health Expectancies ****************/
-void evsij(double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int cij, int estepm,char strstart[] )
+ void evsij(double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int cij, int estepm,char strstart[], int nres )
{
/* Health expectancies, no variances */
@@ -4935,7 +4981,7 @@ void evsij(double ***eij, double x[], in
/* Computed by stepm unit matrices, product of hstepma matrices, stored
in an array of nhstepma length: nhstepma=10, hstepm=4, stepm=6 months */
- hpxij(p3mat,nhstepma,age,hstepm,x,nlstate,stepm,oldm, savm, cij);
+ hpxij(p3mat,nhstepma,age,hstepm,x,nlstate,stepm,oldm, savm, cij, nres);
hf=hstepm*stepm/YEARM; /* Duration of hstepm expressed in year unit. */
@@ -4970,7 +5016,7 @@ void evsij(double ***eij, double x[], in
}
-void cvevsij(double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int cij, int estepm,double delti[],double **matcov,char strstart[] )
+ void cvevsij(double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int cij, int estepm,double delti[],double **matcov,char strstart[], int nres )
{
/* Covariances of health expectancies eij and of total life expectancies according
@@ -5083,8 +5129,8 @@ void cvevsij(double ***eij, double x[],
xp[i] = x[i] + (i==theta ?delti[theta]:0);
xm[i] = x[i] - (i==theta ?delti[theta]:0);
}
- hpxij(p3matp,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, cij);
- hpxij(p3matm,nhstepm,age,hstepm,xm,nlstate,stepm,oldm,savm, cij);
+ hpxij(p3matp,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, cij, nres);
+ hpxij(p3matm,nhstepm,age,hstepm,xm,nlstate,stepm,oldm,savm, cij, nres);
for(j=1; j<= nlstate; j++){
for(i=1; i<=nlstate; i++){
@@ -5125,7 +5171,7 @@ void cvevsij(double ***eij, double x[],
}
/* Computing expectancies */
- hpxij(p3matm,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, cij);
+ hpxij(p3matm,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, cij,nres);
for(i=1; i<=nlstate;i++)
for(j=1; j<=nlstate;j++)
for (h=0, eij[i][j][(int)age]=0; h<=nhstepm-1; h++){
@@ -5180,7 +5226,7 @@ void cvevsij(double ***eij, double x[],
}
/************ Variance ******************/
- void varevsij(char optionfilefiname[], double ***vareij, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int *ncvyearp, int ij, int estepm, int cptcov, int cptcod, int popbased, int mobilav, char strstart[])
+ void varevsij(char optionfilefiname[], double ***vareij, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int *ncvyearp, int ij, int estepm, int cptcov, int cptcod, int popbased, int mobilav, char strstart[], int nres)
{
/* Variance of health expectancies */
/* double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double ** savm,double ftolpl);*/
@@ -5306,7 +5352,7 @@ void cvevsij(double ***eij, double x[],
xp[i] = x[i] + (i==theta ?delti[theta]:0);
}
- prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij, nresult);
if (popbased==1) {
if(mobilav ==0){
@@ -5318,7 +5364,7 @@ void cvevsij(double ***eij, double x[],
}
}
- hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij); /* Returns p3mat[i][j][h] for h=1 to nhstepm */
+ hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij,nres); /* Returns p3mat[i][j][h] for h=1 to nhstepm */
for(j=1; j<= nlstate; j++){
for(h=0; h<=nhstepm; h++){
for(i=1, gp[h][j]=0.;i<=nlstate;i++)
@@ -5338,7 +5384,7 @@ void cvevsij(double ***eij, double x[],
for(i=1; i<=npar; i++) /* Computes gradient x - delta */
xp[i] = x[i] - (i==theta ?delti[theta]:0);
- prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp, ij);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp, ij, nresult);
if (popbased==1) {
if(mobilav ==0){
@@ -5350,7 +5396,7 @@ void cvevsij(double ***eij, double x[],
}
}
- hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij);
+ hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij,nres);
for(j=1; j<= nlstate; j++){ /* Sum of wi * eij = e.j */
for(h=0; h<=nhstepm; h++){
@@ -5415,7 +5461,7 @@ void cvevsij(double ***eij, double x[],
/* end ppptj */
/* x centered again */
- prevalim(prlim,nlstate,x,age,oldm,savm,ftolpl,ncvyearp,ij);
+ prevalim(prlim,nlstate,x,age,oldm,savm,ftolpl,ncvyearp,ij, nresult);
if (popbased==1) {
if(mobilav ==0){
@@ -5431,7 +5477,7 @@ void cvevsij(double ***eij, double x[],
computed over hstepm (estepm) matrices product = hstepm*stepm months)
as a weighted average of prlim.
*/
- hpxij(p3mat,nhstepm,age,hstepm,x,nlstate,stepm,oldm,savm, ij);
+ hpxij(p3mat,nhstepm,age,hstepm,x,nlstate,stepm,oldm,savm, ij, nres);
for(j=nlstate+1;j<=nlstate+ndeath;j++){
for(i=1,gmp[j]=0.;i<= nlstate; i++)
gmp[j] += prlim[i][i]*p3mat[i][j][1];
@@ -5494,7 +5540,7 @@ void cvevsij(double ***eij, double x[],
} /* end varevsij */
/************ Variance of prevlim ******************/
- void varprevlim(char fileres[], double **varpl, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int *ncvyearp, int ij, char strstart[])
+ void varprevlim(char fileres[], double **varpl, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int *ncvyearp, int ij, char strstart[], int nres)
{
/* Variance of prevalence limit for each state ij using current parameters x[] and estimates of neighbourhood give by delti*/
/* double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double **savm,double ftolpl);*/
@@ -5537,9 +5583,9 @@ void cvevsij(double ***eij, double x[],
xp[i] = x[i] + (i==theta ?delti[theta]:0);
}
if((int)age==79 ||(int)age== 80 ||(int)age== 81 )
- prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres);
else
- prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres);
for(i=1;i<=nlstate;i++){
gp[i] = prlim[i][i];
mgp[theta][i] = prlim[i][i];
@@ -5547,9 +5593,9 @@ void cvevsij(double ***eij, double x[],
for(i=1; i<=npar; i++) /* Computes gradient */
xp[i] = x[i] - (i==theta ?delti[theta]:0);
if((int)age==79 ||(int)age== 80 ||(int)age== 81 )
- prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres);
else
- prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres);
for(i=1;i<=nlstate;i++){
gm[i] = prlim[i][i];
mgm[theta][i] = prlim[i][i];
@@ -5957,11 +6003,13 @@ void printinghtml(char fileresu[], char
int popforecast, int prevfcast, int backcast, int estepm , \
double jprev1, double mprev1,double anprev1, double dateprev1, \
double jprev2, double mprev2,double anprev2, double dateprev2){
- int jj1, k1, i1, cpt;
+ int jj1, k1, i1, cpt, k4, nres;
fprintf(fichtm,"");
+ fprintf(fichtm,"", model);
fprintf(fichtm,"- \n");
fprintf(fichtm,"
- - Observed frequency between two states (during the period defined between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf): %s (html file)
\n",
jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,subdirfext3(optionfilefiname,"PHTMFR_",".htm"),subdirfext3(optionfilefiname,"PHTMFR_",".htm"));
@@ -5996,16 +6044,27 @@ void printinghtml(char fileresu[], char
if (cptcovn < 1) {m=1;ncodemax[1]=1;}
jj1=0;
+
+ for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k1=1; k1<=m;k1++){
+ if(TKresult[nres]!= k1)
+ continue;
/* for(i1=1; i1<=ncodemax[k1];i1++){ */
jj1++;
if (cptcovn > 0) {
fprintf(fichtm,"
************ Results for covariates");
for (cpt=1; cpt<=cptcoveff;cpt++){
- fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]);
- printf(" V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]);fflush(stdout);
+ fprintf(fichtm," V%d=%d ",Tvresult[nres][cpt],(int)Tresult[nres][cpt]);
+ printf(" V%d=%d ",Tvresult[nres][cpt],Tresult[nres][cpt]);fflush(stdout);
+ /* fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]); */
+ /* printf(" V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]);fflush(stdout); */
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(fichtm," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);fflush(stdout);
+ }
+
/* if(nqfveff+nqtveff 0) */ /* Test to be done */
fprintf(fichtm," ************\n
");
if(invalidvarcomb[k1]){
@@ -6117,13 +6176,22 @@ See page 'Matrix of variance-covariance
if (cptcovn < 1) {m=1;ncodemax[1]=1;}
jj1=0;
+
+ for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k1=1; k1<=m;k1++){
+ if(TKresult[nres]!= k1)
+ continue;
/* for(i1=1; i1<=ncodemax[k1];i1++){ */
jj1++;
if (cptcovn > 0) {
fprintf(fichtm,"
************ Results for covariates");
for (cpt=1; cpt<=cptcoveff;cpt++) /**< cptcoveff number of variables */
- fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]);
+ fprintf(fichtm," V%d=%d ",Tvresult[nres][cpt],Tresult[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," ************\n
");
if(invalidvarcomb[k1]){
@@ -6153,11 +6221,12 @@ void printinggnuplot(char fileresu[], ch
char dirfileres[132],optfileres[132];
char gplotcondition[132];
- int cpt=0,k1=0,i=0,k=0,j=0,jk=0,k2=0,k3=0,ij=0,l=0;
+ int cpt=0,k1=0,i=0,k=0,j=0,jk=0,k2=0,k3=0,k4=0,ij=0, ijp=0, l=0;
int lv=0, vlv=0, kl=0;
int ng=0;
int vpopbased;
int ioffset; /* variable offset for columns */
+ int nres=0; /* Index of resultline */
/* if((ficgp=fopen(optionfilegnuplot,"a"))==NULL) { */
/* printf("Problem with file %s",optionfilegnuplot); */
@@ -6202,9 +6271,14 @@ void printinggnuplot(char fileresu[], ch
strcpy(optfileres,"vpl");
/* 1eme*/
for (cpt=1; cpt<= nlstate ; cpt ++) { /* For each live state */
- for (k1=1; k1<= m && selected(k1) ; k1 ++) { /* For each valid combination of covariate */
- /* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */
- fprintf(ficgp,"\n# 1st: Period (stable) prevalence with CI: 'VPL_' files ");
+ for (k1=1; k1<= m ; k1 ++) /* For each valid combination of covariate */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ /* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */
+ if(TKresult[nres]!= k1)
+ continue;
+ /* We are interested in selected combination by the resultline */
+ printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);
+ fprintf(ficgp,"\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);
for (k=1; k<=cptcoveff; k++){ /* For each covariate k get corresponding value lv for combination k1 */
lv= decodtabm(k1,k,cptcoveff); /* Should be the value of the covariate corresponding to k1 combination */
/* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */
@@ -6212,21 +6286,27 @@ void printinggnuplot(char fileresu[], ch
/* decodtabm(13,3,4)= 2 because h=13 k= 1 1 (2) 2 */
vlv= nbcode[Tvaraff[k]][lv]; /* vlv is the value of the covariate lv, 0 or 1 */
/* For each combination of covariate k1 (V1=1, V3=0), we printed the current covariate k and its value vlv */
+ printf(" V%d=%d ",Tvaraff[k],vlv);
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
+ printf("\n#\n");
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
continue;
}
-
+
fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"V_"),cpt,k1);
fprintf(ficgp,"\n#set out \"V_%s_%d-%d.svg\" \n",optionfilefiname,cpt,k1);
fprintf(ficgp,"set xlabel \"Age\" \n\
-set ylabel \"Probability\" \n \
-set ter svg size 640, 480\n \
+set ylabel \"Probability\" \n \
+set ter svg size 640, 480\n \
plot [%.f:%.f] \"%s\" every :::%d::%d u 1:2 \"%%lf",ageminpar,fage,subdirf2(fileresu,"VPL_"),k1-1,k1-1);
-
+
for (i=1; i<= nlstate ; i ++) {
if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
else fprintf(ficgp," %%*lf (%%*lf)");
@@ -6273,9 +6353,13 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
fprintf(ficgp,"\nset out \n");
} /* k1 */
} /* cpt */
- /*2 eme*/
- for (k1=1; k1<= m ; k1 ++) {
+
+ /*2 eme*/
+ for (k1=1; k1<= m ; k1 ++)
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
fprintf(ficgp,"\n# 2nd: Total life expectancy with CI: 't' files ");
for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
@@ -6285,6 +6369,10 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
vlv= nbcode[Tvaraff[k]][lv];
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
@@ -6322,15 +6410,18 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
} /* state */
} /* vpopbased */
fprintf(ficgp,"\nset out;set out \"%s_%d.svg\"; replot; set out; \n",subdirf2(optionfilefiname,"E_"),k1); /* Buggy gnuplot */
- } /* k1 */
+ } /* 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(TKresult[nres]!= k)
+ continue;
for (cpt=1; cpt<= nlstate ; cpt ++) {
- fprintf(ficgp,"\n# 3d: Life expectancy with EXP_ files: cov=%d state=%d",k1, cpt);
- for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
+ fprintf(ficgp,"\n# 3d: Life expectancy with EXP_ files: combination=%d state=%d",k1, cpt);
+ for (k=1; k<=cptcoveff; k++){ /* For each covariate dummy combination and each value */
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 */
@@ -6338,6 +6429,10 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
vlv= nbcode[Tvaraff[k]][lv];
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ /* for(k=1; k <= ncovds; k++){ */
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
@@ -6368,7 +6463,10 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
/* 4eme */
/* Survival functions (period) from state i in state j by initial state i */
- for (k1=1; k1<= m ; k1 ++) { /* For each multivariate if any */
+ for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
fprintf(ficgp,"\n#\n#\n# Survival functions in state j : 'LIJ_' files, cov=%d state=%d",k1, cpt);
@@ -6380,6 +6478,9 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
vlv= nbcode[Tvaraff[k]][lv];
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
@@ -6410,7 +6511,10 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
/* 5eme */
/* Survival functions (period) from state i in state j by final state j */
- for (k1=1; k1<= m ; k1 ++) { /* For each covariate if any */
+ for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each inital state */
fprintf(ficgp,"\n#\n#\n# Survival functions in state j and all livestates from state i by final state j: 'lij' files, cov=%d state=%d",k1, cpt);
@@ -6422,6 +6526,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
vlv= nbcode[Tvaraff[k]][lv];
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
@@ -6460,7 +6567,10 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
/* 6eme */
/* CV preval stable (period) for each covariate */
- for (k1=1; k1<= m ; k1 ++) { /* For each covariate combination (1 to m=2**k), if any covariate is present */
+ for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
fprintf(ficgp,"\n#\n#\n#CV preval stable (period): 'pij' files, covariatecombination#=%d state=%d",k1, cpt);
@@ -6472,6 +6582,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
vlv= nbcode[Tvaraff[k]][lv];
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
@@ -6503,7 +6616,10 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
/* 7eme */
if(backcast == 1){
/* CV back preval stable (period) for each covariate */
- for (k1=1; k1<= m ; k1 ++) { /* For each covariate combination (1 to m=2**k), if any covariate is present */
+ for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
fprintf(ficgp,"\n#\n#\n#CV Back preval stable (period): 'pij' files, covariatecombination#=%d state=%d",k1, cpt);
for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
@@ -6514,6 +6630,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
vlv= nbcode[Tvaraff[k]][lv];
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
@@ -6550,7 +6669,10 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
if(prevfcast==1){
/* Projection from cross-sectional to stable (period) for each covariate */
- for (k1=1; k1<= m ; k1 ++) { /* For each covariate combination (1 to m=2**k), if any covariate is present */
+ for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= k1)
+ continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
fprintf(ficgp,"\n#\n#\n#Projection of prevalence to stable (period): 'PROJ_' files, covariatecombination#=%d state=%d",k1, cpt);
for (k=1; k<=cptcoveff; k++){ /* For each correspondig covariate value */
@@ -6561,6 +6683,9 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
vlv= nbcode[Tvaraff[k]][lv];
fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficgp,"\n#\n");
if(invalidvarcomb[k1]){
fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);
@@ -6672,10 +6797,19 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
fprintf(ficgp,"# +exp(a14+b14*age+c14age*age+d14*V1+e14*V1*age)+...)\n");
fprintf(ficgp,"#\n");
for(ng=1; ng<=3;ng++){ /* Number of graphics: first is logit, 2nd is probabilities, third is incidences per year*/
+ fprintf(ficgp,"#Number of graphics: first is logit, 2nd is probabilities, third is incidences per year \n");
+ fprintf(ficgp,"#model=%s \n",model);
fprintf(ficgp,"# ng=%d\n",ng);
- fprintf(ficgp,"# jk=1 to 2^%d=%d\n",cptcoveff,m);
- for(jk=1; jk <=m; jk++) {
- fprintf(ficgp,"# jk=%d\n",jk);
+ fprintf(ficgp,"# jk=1 to 2^%d=%d\n",cptcoveff,m);/* to be checked */
+ for(jk=1; jk <=m; jk++) /* For each combination of covariate */
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ if(TKresult[nres]!= jk)
+ continue;
+ fprintf(ficgp,"# Combination of dummy jk=%d and ",jk);
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
+ fprintf(ficgp,"\n#\n");
fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" ",subdirf2(optionfilefiname,"PE_"),jk,ng);
fprintf(ficgp,"\nset ter svg size 640, 480 ");
if (ng==1){
@@ -6716,18 +6850,49 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
break;
}
ij=1;/* To be checked else nbcode[0][0] wrong */
- for(j=3; j <=ncovmodel-nagesqr; j++) {
+ ijp=1; /* product no age */
+ /* for(j=3; j <=ncovmodel-nagesqr; j++) { */
+ for(j=1; j <=cptcovt; j++) { /* For each covariate of the simplified model */
/* printf("Tage[%d]=%d, j=%d\n", ij, Tage[ij], j); */
- if(ij <=cptcovage) { /* Bug valgrind */
- if((j-2)==Tage[ij]) { /* Bug valgrind */
- fprintf(ficgp,"+p%d*%d*x",i+j+nagesqr-1,nbcode[Tvar[j-2]][codtabm(jk,j-2)]);
- /* fprintf(ficgp,"+p%d*%d*x",i+j+nagesqr-1,nbcode[Tvar[j-2]][codtabm(jk,Tvar[j-2])]); */
+ if(j==Tage[ij]) { /* Product by age */
+ if(ij <=cptcovage) { /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, 2 V5 and V1 */
+ if(Dummy[j]==0){
+ fprintf(ficgp,"+p%d*%d*x",i+j+2+nagesqr-1,Tinvresult[nres][Tvar[j]]);;
+ }else{ /* quantitative */
+ fprintf(ficgp,"+p%d*%f*x",i+j+2+nagesqr-1,Tqinvresult[nres][Tvar[j]]); /* Tqinvresult in decoderesult */
+ /* fprintf(ficgp,"+p%d*%d*x",i+j+nagesqr-1,nbcode[Tvar[j-2]][codtabm(jk,Tvar[j-2])]); */
+ }
ij++;
}
- }
- else
- fprintf(ficgp,"+p%d*%d",i+j+nagesqr-1,nbcode[Tvar[j-2]][codtabm(jk,j-2)]); /* Valgrind bug nbcode */
- }
+ }else if(j==Tprod[ijp]) { /* */
+ /* printf("Tprod[%d]=%d, j=%d\n", ij, Tprod[ijp], j); */
+ if(ijp <=cptcovprod) { /* Product */
+ if(Dummy[Tvard[ijp][1]]==0){/* Vn is dummy */
+ if(Dummy[Tvard[ijp][2]]==0){/* Vn and Vm are dummy */
+ /* fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,nbcode[Tvard[ijp][1]][codtabm(jk,j)],nbcode[Tvard[ijp][2]][codtabm(jk,j)]); */
+ fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tinvresult[nres][Tvard[ijp][2]]);
+ }else{ /* Vn is dummy and Vm is quanti */
+ /* fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,nbcode[Tvard[ijp][1]][codtabm(jk,j)],Tqinvresult[nres][Tvard[ijp][2]]); */
+ fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);
+ }
+ }else{ /* Vn*Vm Vn is quanti */
+ if(Dummy[Tvard[ijp][2]]==0){
+ fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][2]],Tqinvresult[nres][Tvard[ijp][1]]);
+ }else{ /* Both quanti */
+ fprintf(ficgp,"+p%d*%f*%f",i+j+2+nagesqr-1,Tqinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);
+ }
+ }
+ }
+ } else{ /* simple covariate */
+ /* fprintf(ficgp,"+p%d*%d",i+j+2+nagesqr-1,nbcode[Tvar[j]][codtabm(jk,j)]); /\* Valgrind bug nbcode *\/ */
+ if(Dummy[j]==0){
+ fprintf(ficgp,"+p%d*%d",i+j+2+nagesqr-1,Tinvresult[nres][Tvar[j]]); /* */
+ }else{ /* quantitative */
+ fprintf(ficgp,"+p%d*%f",i+j+2+nagesqr-1,Tqinvresult[nres][Tvar[j]]); /* */
+ /* fprintf(ficgp,"+p%d*%d*x",i+j+nagesqr-1,nbcode[Tvar[j-2]][codtabm(jk,Tvar[j-2])]); */
+ }
+ } /* end simple */
+ } /* end j */
}else{
i=i-ncovmodel;
if(ng !=1 ) /* For logit formula of log p11 is more difficult to get */
@@ -6745,8 +6910,8 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
ij=1;
for(j=3; j <=ncovmodel-nagesqr; j++){
- if(ij <=cptcovage) { /* Bug valgrind */
- if((j-2)==Tage[ij]) { /* Bug valgrind */
+ if((j-2)==Tage[ij]) { /* Bug valgrind */
+ if(ij <=cptcovage) { /* Bug valgrind */
fprintf(ficgp,"+p%d*%d*x",k3+(k1-1)*ncovmodel+1+j-2+nagesqr,nbcode[Tvar[j-2]][codtabm(jk,j-2)]);
/* fprintf(ficgp,"+p%d*%d*x",k3+(k1-1)*ncovmodel+1+j-2+nagesqr,nbcode[Tvar[j-2]][codtabm(jk,Tvar[j-2])]); */
ij++;
@@ -6898,7 +7063,7 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
} /* end bad */
for (age=bage; age<=fage; age++){
- printf("%d %d ", cptcod, (int)age);
+ /* printf("%d %d ", cptcod, (int)age); */
sumnewp[cptcod]=0.;
sumnewm[cptcod]=0.;
for (i=1; i<=nlstate;i++){
@@ -6937,13 +7102,13 @@ plot [%.f:%.f] ", ageminpar, agemaxpar)
/************** Forecasting ******************/
-void prevforecast(char fileres[], double anproj1, double mproj1, double jproj1, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double bage, double fage, int firstpass, int lastpass, double anproj2, double p[], int cptcoveff){
+ void prevforecast(char fileres[], double anproj1, double mproj1, double jproj1, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double bage, double fage, int firstpass, int lastpass, double anproj2, double p[], int cptcoveff){
/* proj1, year, month, day of starting projection
agemin, agemax range of age
dateprev1 dateprev2 range of dates during which prevalence is computed
anproj2 year of en of projection (same day and month as proj1).
*/
- int yearp, stepsize, hstepm, nhstepm, j, k, cptcod, i, h, i1;
+ int yearp, stepsize, hstepm, nhstepm, j, k, cptcod, i, h, i1, k4, nres=0;
double agec; /* generic age */
double agelim, ppij, yp,yp1,yp2,jprojmean,mprojmean,anprojmean;
double *popeffectif,*popcount;
@@ -6965,8 +7130,8 @@ void prevforecast(char fileres[], double
printf("Problem with forecast resultfile: %s\n", fileresf);
fprintf(ficlog,"Problem with forecast resultfile: %s\n", fileresf);
}
- printf("Computing forecasting: result on file '%s', please wait... \n", fileresf);
- fprintf(ficlog,"Computing forecasting: result on file '%s', please wait... \n", fileresf);
+ printf("\nComputing forecasting: result on file '%s', please wait... \n", fileresf);
+ fprintf(ficlog,"\nComputing forecasting: result on file '%s', please wait... \n", fileresf);
if (cptcoveff==0) ncodemax[cptcoveff]=1;
@@ -6997,7 +7162,10 @@ void prevforecast(char fileres[], double
fprintf(ficresf,"#****** Routine prevforecast **\n");
/* if (h==(int)(YEARM*yearp)){ */
- for(k=1;k<=i1;k++){
+ for(nres=1; nres <= nresult; nres++) /* For each resultline */
+ for(k=1; k<=i1;k++){
+ if(TKresult[nres]!= k)
+ continue;
if(invalidvarcomb[k]){
printf("\nCombination (%d) projection ignored because no cases \n",k);
continue;
@@ -7006,6 +7174,10 @@ void prevforecast(char fileres[], double
for(j=1;j<=cptcoveff;j++) {
fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficresf," yearproj age");
for(j=1; j<=nlstate+ndeath;j++){
for(i=1; i<=nlstate;i++)
@@ -7020,7 +7192,7 @@ void prevforecast(char fileres[], double
nhstepm = nhstepm/hstepm;
p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
oldm=oldms;savm=savms;
- hpxij(p3mat,nhstepm,agec,hstepm,p,nlstate,stepm,oldm,savm, k);
+ hpxij(p3mat,nhstepm,agec,hstepm,p,nlstate,stepm,oldm,savm, k,nres);
for (h=0; h<=nhstepm; h++){
if (h*hstepm/YEARM*stepm ==yearp) {
@@ -7890,7 +8062,7 @@ int readdata(char datafile[], int firsto
void removefirstspace(char **stri){/*, char stro[]) {*/
char *p1 = *stri, *p2 = *stri;
- if (*p2 == ' ')
+ while (*p2 == ' ')
p2++;
/* while ((*p1++ = *p2++) !=0) */
/* ; */
@@ -7901,10 +8073,10 @@ void removefirstspace(char **stri){/*, c
*stri=p2;
}
-int decoderesult ( char resultline[])
+int decoderesult ( char resultline[], int nres)
/**< This routine decode one result line and returns the combination # of dummy covariates only **/
{
- int j=0, k=0, k1=0, k2=0, match=0;
+ 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];
@@ -7942,13 +8114,13 @@ int decoderesult ( char resultline[])
if (nbocc(stra,'=') >0)
strcpy(resultsav,stra); /* and analyzes it */
}
- /* Checking if no missing or useless values in comparison of current model needs */
- for(k1=1; k1<= cptcovt ;k1++){ /* model line */
- if(Typevar[k1]==0){
+ /* Checking for missing or useless values in comparison of current model needs */
+ for(k1=1; k1<= cptcovt ;k1++){ /* model line V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
+ if(Typevar[k1]==0){ /* Single covariate in model */
match=0;
- for(k2=1; k2 <=j;k2++){
- if(Tvar[k1]==Tvarsel[k2]) {
- modelresult[k2]=k1;
+ for(k2=1; k2 <=j;k2++){/* result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */
+ if(Tvar[k1]==Tvarsel[k2]) {/* Tvar[1]=5 == Tvarsel[2]=5 */
+ modelresult[k2]=k1;/* modelresult[2]=1 modelresult[1]=2 modelresult[3]=3 modelresult[6]=4 modelresult[9]=5 */
match=1;
break;
}
@@ -7958,13 +8130,13 @@ int decoderesult ( char resultline[])
}
}
}
-
- for(k2=1; k2 <=j;k2++){ /* result line */
+ /* Checking for missing or useless values in comparison of current model needs */
+ for(k2=1; k2 <=j;k2++){ /* result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */
match=0;
- for(k1=1; k1<= cptcovt ;k1++){ /* model line */
- if(Typevar[k1]==0){
- if(Tvar[k1]==Tvarsel[k2]) {
- resultmodel[k1]=k2;
+ for(k1=1; k1<= cptcovt ;k1++){ /* model line V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
+ if(Typevar[k1]==0){ /* Single */
+ if(Tvar[k1]==Tvarsel[k2]) { /* Tvar[2]=4 == Tvarsel[1]=4 */
+ resultmodel[k1]=k2; /* resultmodel[2]=1 resultmodel[1]=2 resultmodel[3]=3 resultmodel[6]=4 resultmodel[9]=5 */
++match;
}
}
@@ -7975,17 +8147,50 @@ int decoderesult ( char resultline[])
printf("Error in result line: %d doubled; result: %s, model=%s\n",k2, resultline, model);
}
}
-
+
/* 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 */
+ /* result line V4=1 V5=25.1 V3=0 V2=8 V1=1 */
+ /* should give a combination of dummy V4=1, V3=0, V1=1 => V4*2**(0) + V3*2**(1) + V1*2**(2) = 5 + (1offset) = 6*/
+ /* result line V4=1 V5=24.1 V3=1 V2=8 V1=0 */
+ /* should give a combination of dummy V4=1, V3=1, V1=0 => V4*2**(0) + V3*2**(1) + V1*2**(2) = 3 + (1offset) = 4*/
+ /* 1 0 0 0 */
+ /* 2 1 0 0 */
+ /* 3 0 1 0 */
+ /* 4 1 1 0 */ /* V4=1, V3=1, V1=0 */
+ /* 5 0 0 1 */
+ /* 6 1 0 1 */ /* V4=1, V3=0, V1=1 */
+ /* 7 0 1 1 */
+ /* 8 1 1 1 */
+ /* V(Tvresult)=Tresult V4=1 V3=0 V1=1 Tresult[nres=1][2]=0 */
+ /* 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++){ /* model line */
+ if( Dummy[k1]==0 && Typevar[k1]==0 ){ /* Single dummy */
+ k3= resultmodel[k1]; /* resultmodel[2(V4)] = 1=k3 */
+ k2=(int)Tvarsel[k3]; /* Tvarsel[resultmodel[2]]= Tvarsel[1] = 4=k2 */
+ k+=Tvalsel[k3]*pow(2,k4); /* Tvalsel[1]=1 */
+ Tresult[nres][k4+1]=Tvalsel[k3];/* Tresult[nres][1]=1(V4=1) Tresult[nres][2]=0(V3=0) */
+ Tvresult[nres][k4+1]=(int)Tvarsel[k3];/* Tvresult[nres][1]=4 Tvresult[nres][3]=1 */
+ Tinvresult[nres][(int)Tvarsel[k3]]=Tvalsel[k3]; /* Tinvresult[nres][4]=1 */
+ printf("Decoderesult Dummy k=%d, V(k2=V%d)= Tvalsel[%d]=%d, 2**(%d)\n",k, k2, k3, (int)Tvalsel[k3], k4);
+ k4++;;
+ } else if( Dummy[k1]==1 && Typevar[k1]==0 ){ /* Single quantitative */
+ k3q= resultmodel[k1]; /* resultmodel[2] = 1=k3 */
+ k2q=(int)Tvarsel[k3q]; /* Tvarsel[resultmodel[2]]= 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 */
+ Tqinvresult[nres][(int)Tvarsel[k3q]]=Tvalsel[k3q]; /* Tqinvresult[nres][5]=25.1 */
+ printf("Decoderesult Quantitative nres=%d, V(k2q=V%d)= Tvalsel[%d]=%d, Tvarsel[%d]=%f\n",nres, k2q, k3q, Tvarsel[k3q], k3q, Tvalsel[k3q]);
+ k4q++;;
+ }
+ }
+ TKresult[nres]=++k; /* Combination for the nresult and the model */
return (0);
}
-int selected( int kvar){ /* Selected combination of covariates */
- if(Tvarsel[kvar])
- return (0);
- else
- return(1);
-}
+
int decodemodel( char model[], int lastobs)
/**< This routine decodes the model and returns:
* Model V1+V2+V3+V8+V7*V8+V5*V6+V8*age+V3*age+age*age
@@ -8795,7 +9000,7 @@ void syscompilerinfo(int logged)
int prevalence_limit(double *p, double **prlim, double ageminpar, double agemaxpar, double ftolpl, int *ncvyearp){
/*--------------- Prevalence limit (period or stable prevalence) --------------*/
- int i, j, k, i1 ;
+ int i, j, k, i1, k4=0, nres=0 ;
/* double ftolpl = 1.e-10; */
double age, agebase, agelim;
double tot;
@@ -8823,7 +9028,11 @@ int prevalence_limit(double *p, double *
i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */
if (cptcovn < 1){i1=1;}
+ for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k=1; k<=i1;k++){
+ if(TKresult[nres]!= k)
+ continue;
+
/* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
/* for(cptcov=1,k=0;cptcov<=1;cptcov++){ */
//for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){
@@ -8838,6 +9047,10 @@ int prevalence_limit(double *p, double *
printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
}
+ for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ }
fprintf(ficrespl,"******\n");
printf("******\n");
fprintf(ficlog,"******\n");
@@ -8857,7 +9070,7 @@ int prevalence_limit(double *p, double *
for (age=agebase; age<=agelim; age++){
/* for (age=agebase; age<=agebase; age++){ */
- prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, ncvyearp, k);
+ prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, ncvyearp, k, nres);
fprintf(ficrespl,"%.0f ",age );
for(j=1;j<=cptcoveff;j++)
fprintf(ficrespl,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
@@ -8879,7 +9092,7 @@ int back_prevalence_limit(double *p, dou
/* Computes the back prevalence limit for any combination of covariate values
* at any age between ageminpar and agemaxpar
*/
- int i, j, k, i1 ;
+ int i, j, k, i1, nres=0 ;
/* double ftolpl = 1.e-10; */
double age, agebase, agelim;
double tot;
@@ -8910,7 +9123,10 @@ int back_prevalence_limit(double *p, dou
i1=pow(2,cptcoveff);
if (cptcovn < 1){i1=1;}
- for(k=1; k<=i1;k++){
+ for(nres=1; nres <= nresult; nres++) /* For each resultline */
+ for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
+ if(TKresult[nres]!= k)
+ continue;
//printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
fprintf(ficresplb,"#******");
printf("#******");
@@ -8920,6 +9136,11 @@ int back_prevalence_limit(double *p, dou
printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
}
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficresplb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
fprintf(ficresplb,"******\n");
printf("******\n");
fprintf(ficlog,"******\n");
@@ -8978,7 +9199,7 @@ int hPijx(double *p, int bage, int fage)
int agelim;
int hstepm;
int nhstepm;
- int h, i, i1, j, k;
+ int h, i, i1, j, k, k4, nres=0;
double agedeb;
double ***p3mat;
@@ -9005,10 +9226,17 @@ int hPijx(double *p, int bage, int fage)
/* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
/* /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */
/* k=k+1; */
- for (k=1; k <= (int) pow(2,cptcoveff); k++){
+ for(nres=1; nres <= nresult; nres++) /* For each resultline */
+ for(k=1; k<=i1;k++){
+ if(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,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 */
@@ -9019,7 +9247,7 @@ int hPijx(double *p, int bage, int fage)
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);
+ 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++)
@@ -9145,6 +9373,7 @@ int main(int argc, char *argv[])
int itimes;
int NDIM=2;
int vpopbased=0;
+ int nres=0;
char ca[32], cb[32];
/* FILE *fichtm; *//* Html File */
@@ -10548,6 +10777,7 @@ Please run with mle=-1 to get a correct
/* day and month of proj2 are not used but only year anproj2.*/
/* Results */
+ nresult=0;
while(fgets(line, MAXLINE, ficpar)) {
/* If line starts with a # it is a comment */
if (line[0] == '#') {
@@ -10565,8 +10795,14 @@ Please run with mle=-1 to get a correct
else if (num_filled != 1){
printf("ERROR %d: result line should be at minimum 'result=' %s\n",num_filled, line);
}
- printf("Result %d: result line should be at minimum 'line=' %s, result=%s\n",num_filled, line, resultline);
- decoderesult(resultline);
+ nresult++; /* Sum of resultlines */
+ printf("Result %d: result=%s\n",nresult, resultline);
+ if(nresult > MAXRESULTLINES){
+ printf("ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
+ fprintf(ficlog,"ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
+ goto end;
+ }
+ decoderesult(resultline, nresult); /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */
while(fgets(line, MAXLINE, ficpar)) {
/* If line starts with a # it is a comment */
if (line[0] == '#') {
@@ -10653,6 +10889,7 @@ Please run with mle=-1 to get a correct
mobaverages[i][j][k]=0.;
mobaverage=mobaverages;
if (mobilav!=0) {
+ printf("Movingaveraging observed prevalence\n");
if (movingaverage(probs, ageminpar, agemaxpar, mobaverage, mobilav)!=0){
fprintf(ficlog," Error in movingaverage mobilav=%d\n",mobilav);
printf(" Error in movingaverage mobilav=%d\n",mobilav);
@@ -10661,6 +10898,7 @@ Please run with mle=-1 to get a correct
/* /\* Prevalence for each covariates in probs[age][status][cov] *\/ */
/* prevalence(probs, ageminpar, agemaxpar, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass); */
else if (mobilavproj !=0) {
+ printf("Movingaveraging projected observed prevalence\n");
if (movingaverage(probs, ageminpar, agemaxpar, mobaverage, mobilavproj)!=0){
fprintf(ficlog," Error in movingaverage mobilavproj=%d\n",mobilavproj);
printf(" Error in movingaverage mobilavproj=%d\n",mobilavproj);
@@ -10716,16 +10954,29 @@ Please run with mle=-1 to get a correct
printf("Computing Health Expectancies: result on file '%s' ...", filerese);fflush(stdout);
fprintf(ficlog,"Computing Health Expectancies: result on file '%s' ...", filerese);fflush(ficlog);
- for (k=1; k <= (int) pow(2,cptcoveff); k++){ /* For any combination of dummy covariates, fixed and varying */
+ 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(TKresult[nres]!= k)
+ continue;
fprintf(ficreseij,"\n#****** ");
+ printf("\n#****** ");
for(j=1;j<=cptcoveff;j++) {
fprintf(ficreseij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ }
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficreseij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
}
fprintf(ficreseij,"******\n");
+ printf("******\n");
eij=ma3x(1,nlstate,1,nlstate,(int) bage, (int) fage);
oldm=oldms;savm=savms;
- evsij(eij, p, nlstate, stepm, (int) bage, (int)fage, oldm, savm, k, estepm, strstart);
+ evsij(eij, p, nlstate, stepm, (int) bage, (int)fage, oldm, savm, k, estepm, strstart, nres);
free_ma3x(eij,1,nlstate,1,nlstate,(int) bage, (int)fage);
}
@@ -10776,15 +11027,26 @@ Please run with mle=-1 to get a correct
/*for(cptcov=1,k=0;cptcov<=i1;cptcov++){
for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*/
- for (k=1; k <= (int) pow(2,cptcoveff); k++){
- printf("\n#****** ");
- fprintf(ficrest,"\n#****** ");
- fprintf(ficlog,"\n#****** ");
+ 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(TKresult[nres]!= k)
+ continue;
+ printf("\n#****** Selected:");
+ fprintf(ficrest,"\n#****** Selected:");
+ fprintf(ficlog,"\n#****** Selected:");
for(j=1;j<=cptcoveff;j++){
printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
fprintf(ficrest,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
fprintf(ficlog,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
}
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficrest," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
fprintf(ficrest,"******\n");
fprintf(ficlog,"******\n");
printf("******\n");
@@ -10795,19 +11057,26 @@ Please run with mle=-1 to get a correct
fprintf(ficresstdeij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
fprintf(ficrescveij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
}
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ fprintf(ficresstdeij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficrescveij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
fprintf(ficresstdeij,"******\n");
fprintf(ficrescveij,"******\n");
fprintf(ficresvij,"\n#****** ");
for(j=1;j<=cptcoveff;j++)
fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ fprintf(ficresvij," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
fprintf(ficresvij,"******\n");
eij=ma3x(1,nlstate,1,nlstate,(int) bage, (int) fage);
oldm=oldms;savm=savms;
- printf(" cvevsij combination#=%d, ",k);
- fprintf(ficlog, " cvevsij combination#=%d, ",k);
- cvevsij(eij, p, nlstate, stepm, (int) bage, (int)fage, oldm, savm, k, estepm, delti, matcov, strstart);
+ printf(" cvevsij ");
+ fprintf(ficlog, " cvevsij ");
+ cvevsij(eij, p, nlstate, stepm, (int) bage, (int)fage, oldm, savm, k, estepm, delti, matcov, strstart, nres);
printf(" end cvevsij \n ");
fprintf(ficlog, " end cvevsij \n ");
@@ -10824,7 +11093,7 @@ Please run with mle=-1 to get a correct
cptcod= 0; /* To be deleted */
printf("varevsij vpopbased=%d \n",vpopbased);
fprintf(ficlog, "varevsij vpopbased=%d \n",vpopbased);
- varevsij(optionfilefiname, vareij, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl, &ncvyear, k, estepm, cptcov,cptcod,vpopbased,mobilav, strstart); /* cptcod not initialized Intel */
+ varevsij(optionfilefiname, vareij, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl, &ncvyear, k, estepm, cptcov,cptcod,vpopbased,mobilav, strstart, nres); /* cptcod not initialized Intel */
fprintf(ficrest,"# Total life expectancy with std error and decomposition into time to be expected in each health state\n# (weighted average of eij where weights are ");
if(vpopbased==1)
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);
@@ -10838,7 +11107,7 @@ Please run with mle=-1 to get a correct
printf("Computing age specific period (stable) prevalences in each health state \n");
fprintf(ficlog,"Computing age specific period (stable) prevalences in each health state \n");
for(age=bage; age <=fage ;age++){
- prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, &ncvyear, k); /*ZZ Is it the correct prevalim */
+ prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, &ncvyear, k, nres); /*ZZ Is it the correct prevalim */
if (vpopbased==1) {
if(mobilav ==0){
for(i=1; i<=nlstate;i++)
@@ -10875,11 +11144,11 @@ Please run with mle=-1 to get a correct
free_ma3x(eij,1,nlstate,1,nlstate,(int) bage, (int)fage);
free_ma3x(vareij,1,nlstate,1,nlstate,(int) bage, (int)fage);
free_vector(epj,1,nlstate+1);
- printf("done \n");fflush(stdout);
- fprintf(ficlog,"done\n");fflush(ficlog);
+ printf("done selection\n");fflush(stdout);
+ fprintf(ficlog,"done selection\n");fflush(ficlog);
/*}*/
- } /* End k */
+ } /* End k selection */
printf("done State-specific expectancies\n");fflush(stdout);
fprintf(ficlog,"done State-specific expectancies\n");fflush(ficlog);
@@ -10898,7 +11167,13 @@ Please run with mle=-1 to get a correct
/*for(cptcov=1,k=0;cptcov<=i1;cptcov++){
for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*/
- for (k=1; k <= (int) pow(2,cptcoveff); k++){
+ 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(TKresult[nres]!= k)
+ continue;
fprintf(ficresvpl,"\n#****** ");
printf("\n#****** ");
fprintf(ficlog,"\n#****** ");
@@ -10907,13 +11182,18 @@ Please run with mle=-1 to get a correct
fprintf(ficlog,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
}
+ for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
+ printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficresvpl," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
+ }
fprintf(ficresvpl,"******\n");
printf("******\n");
fprintf(ficlog,"******\n");
varpl=matrix(1,nlstate,(int) bage, (int) fage);
oldm=oldms;savm=savms;
- varprevlim(fileres, varpl, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl, &ncvyear, k, strstart);
+ varprevlim(fileres, varpl, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl, &ncvyear, k, strstart, nres);
free_matrix(varpl,1,nlstate,(int) bage, (int)fage);
/*}*/
}