/* $Id$
$State$
$Log$
+ Revision 1.241 2016/08/29 17:17:25 brouard
+ Summary: gnuplot problem in Back projection to fix
+
Revision 1.240 2016/08/29 07:53:18 brouard
Summary: Better
/* double **bprevalim(double **bprlim, double ***prevacurrent, int nlstate, double x[], double age, double ageminpar, double agemaxpar, double **oldm, double **savm, double **dnewm, double **doldm, double **dsavm, double ftolpl, int *ncvyear, int ij) */
/* double **bprevalim(double **bprlim, double ***prevacurrent, int nlstate, double x[], double age, double **oldm, double **savm, double **dnewm, double **doldm, double **dsavm, double ftolpl, int *ncvyear, int ij) */
- double **bprevalim(double **bprlim, double ***prevacurrent, int nlstate, double x[], double age, double ftolpl, int *ncvyear, int ij)
+ double **bprevalim(double **bprlim, double ***prevacurrent, int nlstate, double x[], double age, double ftolpl, int *ncvyear, int ij, int nres)
{
/* Computes the prevalence limit in each live state at age x and covariate ij by left multiplying the unit
matrix by transitions matrix until convergence is reached with precision ftolpl */
cov[2]=agefin;
if(nagesqr==1)
cov[3]= agefin*agefin;;
- for (k=1; k<=cptcovn;k++) {
- /* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,Tvar[k])]; */
- 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])]); */
+ 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("bprevalim 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<=cptcovn;k++) { */
+ /* /\* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,Tvar[k])]; *\/ */
+ /* 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])]); *\/ */
+ /* } */
+ 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("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]); */
+ }
+ /* 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])]; *\/ */
+ /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)] * nbcode[Tvard[k][2]][codtabm(ij,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]]=Tqresult[nres][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++){ /* 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]];
+ }
+ }
}
- 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])]; */
- cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)] * nbcode[Tvard[k][2]][codtabm(ij,k)];
/*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]);*/
/* }else */
doldm[ii][j]=(ii==j ? 1./sumnew : 0.0);
}else{
- printf("ii=%d, i=%d, doldm=%lf dsavm=%lf, probs=%lf, sumnew=%lf,agefin=%d\n",ii,j,doldm[ii][j],dsavm[ii][j],prevacurrent[(int)agefin][ii][ij],sumnew, (int)agefin);
+ ;
+ /* printf("ii=%d, i=%d, doldm=%lf dsavm=%lf, probs=%lf, sumnew=%lf,agefin=%d\n",ii,j,doldm[ii][j],dsavm[ii][j],prevacurrent[(int)agefin][ii][ij],sumnew, (int)agefin); */
}
} /*End ii */
} /* End j, At the end doldm is diag[1/(w_1p1i+w_2 p2i)] */
*/
for(mi=1; mi<= wav[i]-1; mi++){
for(k=1; k <= ncovv ; k++){ /* Varying covariates (single and product but no age )*/
- cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]][i];
+ /* cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]][i]; */
+ cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]-ncovcol-nqv][i];
}
for (ii=1;ii<=nlstate+ndeath;ii++)
for (j=1;j<=nlstate+ndeath;j++){
if(nagesqr==1)
cov[3]= agexact*agexact; /* Should be changed here */
for (kk=1; kk<=cptcovage;kk++) {
+ if(!FixedV[Tvar[Tage[kk]]])
cov[Tage[kk]+2+nagesqr]=covar[Tvar[Tage[kk]]][i]*agexact; /* Tage[kk] gives the data-covariate associated with age */
+ else
+ cov[Tage[kk]+2+nagesqr]=cotvar[mw[mi][i]][Tvar[Tage[kk]]-ncovcol-nqv][i]*agexact;
}
out=matprod2(newm,oldm,1,nlstate+ndeath,1,nlstate+ndeath,
1,nlstate+ndeath,pmij(pmmij,cov,ncovmodel,x,nlstate));
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 (single and product but no age )*/
- cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]][i];
- }
+ /* cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]][i]; */
+ cov[ioffset+TvarVind[k]]=cotvar[mw[mi][i]][Tvar[TvarVind[k]]-ncovcol-nqv][i];
+ }
/* for(itv=1; itv <= ntveff; itv++){ /\* Varying dummy covariates (single??)*\/ */
- /* iv= Tvar[Tmodelind[ioffset-2-nagesqr-cptcovage+itv]]-ncovcol-nqv; /\* Counting the # varying covariate from 1 to ntveff *\/ */
- /* cov[ioffset+iv]=cotvar[mw[mi][i]][iv][i]; */
- /* k=ioffset-2-nagesqr-cptcovage+itv; /\* position in simple model *\/ */
- /* cov[ioffset+itv]=cotvar[mw[mi][i]][TmodelInvind[itv]][i]; */
- /* printf(" i=%d,mi=%d,itv=%d,TmodelInvind[itv]=%d,cotvar[mw[mi][i]][TmodelInvind[itv]][i]=%f\n", i, mi, itv, TmodelInvind[itv],cotvar[mw[mi][i]][TmodelInvind[itv]][i]); */
+ /* iv= Tvar[Tmodelind[ioffset-2-nagesqr-cptcovage+itv]]-ncovcol-nqv; /\* Counting the # varying covariate from 1 to ntveff *\/ */
+ /* cov[ioffset+iv]=cotvar[mw[mi][i]][iv][i]; */
+ /* k=ioffset-2-nagesqr-cptcovage+itv; /\* position in simple model *\/ */
+ /* cov[ioffset+itv]=cotvar[mw[mi][i]][TmodelInvind[itv]][i]; */
+ /* printf(" i=%d,mi=%d,itv=%d,TmodelInvind[itv]=%d,cotvar[mw[mi][i]][TmodelInvind[itv]][i]=%f\n", i, mi, itv, TmodelInvind[itv],cotvar[mw[mi][i]][TmodelInvind[itv]][i]); */
/* for(iqtv=1; iqtv <= nqtveff; iqtv++){ /\* Varying quantitatives covariates *\/ */
- /* iv=TmodelInvQind[iqtv]; /\* Counting the # varying covariate from 1 to ntveff *\/ */
- /* /\* printf(" i=%d,mi=%d,iqtv=%d,TmodelInvQind[iqtv]=%d,cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i]=%f\n", i, mi, iqtv, TmodelInvQind[iqtv],cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i]); *\/ */
- /* cov[ioffset+ntveff+iqtv]=cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i]; */
+ /* iv=TmodelInvQind[iqtv]; /\* Counting the # varying covariate from 1 to ntveff *\/ */
+ /* /\* printf(" i=%d,mi=%d,iqtv=%d,TmodelInvQind[iqtv]=%d,cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i]=%f\n", i, mi, iqtv, TmodelInvQind[iqtv],cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i]); *\/ */
+ /* cov[ioffset+ntveff+iqtv]=cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i]; */
/* } */
for (ii=1;ii<=nlstate+ndeath;ii++)
- for (j=1;j<=nlstate+ndeath;j++){
- oldm[ii][j]=(ii==j ? 1.0 : 0.0);
- savm[ii][j]=(ii==j ? 1.0 : 0.0);
- }
+ for (j=1;j<=nlstate+ndeath;j++){
+ 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<dh[mi][i]; d++){ /* Delay between two effective waves */
- /*dh[m][i] or dh[mw[mi][i]][i] is the delay between two effective waves m=mw[mi][i]
- and mw[mi+1][i]. dh depends on stepm.*/
- newm=savm;
- agexact=agev[mw[mi][i]][i]+d*stepm/YEARM;
- cov[2]=agexact;
- if(nagesqr==1)
- cov[3]= agexact*agexact;
- for (kk=1; kk<=cptcovage;kk++) {
- cov[Tage[kk]+2+nagesqr]=covar[Tvar[Tage[kk]]][i]*agexact;
- }
- /* printf("i=%d,mi=%d,d=%d,mw[mi][i]=%d\n",i, mi,d,mw[mi][i]); */
- /* savm=pmij(pmmij,cov,ncovmodel,x,nlstate); */
- out=matprod2(newm,oldm,1,nlstate+ndeath,1,nlstate+ndeath,
- 1,nlstate+ndeath,pmij(pmmij,cov,ncovmodel,x,nlstate));
- /* out=matprod2(newm,oldm,1,nlstate+ndeath,1,nlstate+ndeath, */
- /* 1,nlstate+ndeath,pmij(pmmij,cov,ncovmodel,x,nlstate)); */
- savm=oldm;
- oldm=newm;
+ /*dh[m][i] or dh[mw[mi][i]][i] is the delay between two effective waves m=mw[mi][i]
+ and mw[mi+1][i]. dh depends on stepm.*/
+ newm=savm;
+ agexact=agev[mw[mi][i]][i]+d*stepm/YEARM;
+ cov[2]=agexact;
+ if(nagesqr==1)
+ cov[3]= agexact*agexact;
+ for (kk=1; kk<=cptcovage;kk++) {
+ if(!FixedV[Tvar[Tage[kk]]])
+ cov[Tage[kk]+2+nagesqr]=covar[Tvar[Tage[kk]]][i]*agexact;
+ else
+ cov[Tage[kk]+2+nagesqr]=cotvar[mw[mi][i]][Tvar[Tage[kk]]-ncovcol-nqv][i]*agexact;
+ }
+ /* printf("i=%d,mi=%d,d=%d,mw[mi][i]=%d\n",i, mi,d,mw[mi][i]); */
+ /* savm=pmij(pmmij,cov,ncovmodel,x,nlstate); */
+ out=matprod2(newm,oldm,1,nlstate+ndeath,1,nlstate+ndeath,
+ 1,nlstate+ndeath,pmij(pmmij,cov,ncovmodel,x,nlstate));
+ /* out=matprod2(newm,oldm,1,nlstate+ndeath,1,nlstate+ndeath, */
+ /* 1,nlstate+ndeath,pmij(pmmij,cov,ncovmodel,x,nlstate)); */
+ savm=oldm;
+ oldm=newm;
} /* end mult */
s1=s[mw[mi][i]][i];
* is higher than the multiple of stepm and negative otherwise.
*/
if( s2 > nlstate && (mle <5) ){ /* Jackson */
- lli=log(out[s1][s2] - savm[s1][s2]);
+ lli=log(out[s1][s2] - savm[s1][s2]);
} else if ( s2==-1 ) { /* alive */
- for (j=1,survp=0. ; j<=nlstate; j++)
- survp += (1.+bbh)*out[s1][j]- bbh*savm[s1][j];
- lli= log(survp);
+ for (j=1,survp=0. ; j<=nlstate; j++)
+ survp += (1.+bbh)*out[s1][j]- bbh*savm[s1][j];
+ lli= log(survp);
}else if (mle==1){
- lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */
+ lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */
} else if(mle==2){
- lli= (savm[s1][s2]>(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]):log((1.+bbh)*out[s1][s2])); /* linear interpolation */
+ lli= (savm[s1][s2]>(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]):log((1.+bbh)*out[s1][s2])); /* linear interpolation */
} else if(mle==3){ /* exponential inter-extrapolation */
- lli= (savm[s1][s2]>(double)1.e-8 ?(1.+bbh)*log(out[s1][s2])- bbh*log(savm[s1][s2]):log((1.+bbh)*out[s1][s2])); /* exponential inter-extrapolation */
+ lli= (savm[s1][s2]>(double)1.e-8 ?(1.+bbh)*log(out[s1][s2])- bbh*log(savm[s1][s2]):log((1.+bbh)*out[s1][s2])); /* exponential inter-extrapolation */
} else if (mle==4){ /* mle=4 no inter-extrapolation */
- lli=log(out[s1][s2]); /* Original formula */
+ lli=log(out[s1][s2]); /* Original formula */
} else{ /* mle=0 back to 1 */
- lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */
- /*lli=log(out[s1][s2]); */ /* Original formula */
+ lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */
+ /*lli=log(out[s1][s2]); */ /* Original formula */
} /* End of if */
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]); */
if(globpr){
- fprintf(ficresilk,"%9ld %6.1f %6.1f %6d %2d %2d %2d %2d %3d %15.6f %8.4f %8.3f\
+ fprintf(ficresilk,"%9ld %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,out[s1][s2],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);
- }
- fprintf(ficresilk," %10.6f\n", -llt);
+ 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,out[s1][s2],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);
+ }
+ fprintf(ficresilk," %10.6f\n", -llt);
}
} /* end of wave */
} /* end of individual */
/*********** Tricode ****************************/
void tricode(int *cptcov, int *Tvar, int **nbcode, int imx, int *Ndum)
-{
- /**< Uses cptcovn+2*cptcovprod as the number of covariates */
- /* Tvar[i]=atoi(stre); find 'n' in Vn and stores in Tvar. If model=V2+V1 Tvar[1]=2 and Tvar[2]=1
- * Boring subroutine which should only output nbcode[Tvar[j]][k]
- * Tvar[5] in V2+V1+V3*age+V2*V4 is 4 (V4) even it is a time varying or quantitative variable
- * nbcode[Tvar[5]][1]= nbcode[4][1]=0, nbcode[4][2]=1 (usually);
- */
+ {
+ /**< Uses cptcovn+2*cptcovprod as the number of covariates */
+ /* Tvar[i]=atoi(stre); find 'n' in Vn and stores in Tvar. If model=V2+V1 Tvar[1]=2 and Tvar[2]=1
+ * Boring subroutine which should only output nbcode[Tvar[j]][k]
+ * Tvar[5] in V2+V1+V3*age+V2*V4 is 4 (V4) even it is a time varying or quantitative variable
+ * nbcode[Tvar[5]][1]= nbcode[4][1]=0, nbcode[4][2]=1 (usually);
+ */
- int ij=1, k=0, j=0, i=0, maxncov=NCOVMAX;
- int modmaxcovj=0; /* Modality max of covariates j */
- int cptcode=0; /* Modality max of covariates j */
- int modmincovj=0; /* Modality min of covariates j */
+ int ij=1, k=0, j=0, i=0, maxncov=NCOVMAX;
+ int modmaxcovj=0; /* Modality max of covariates j */
+ int cptcode=0; /* Modality max of covariates j */
+ int modmincovj=0; /* Modality min of covariates j */
- /* cptcoveff=0; */
- /* *cptcov=0; */
+ /* cptcoveff=0; */
+ /* *cptcov=0; */
- for (k=1; k <= maxncov; k++) ncodemax[k]=0; /* Horrible constant again replaced by NCOVMAX */
-
- /* Loop on covariates without age and products and no quantitative variable */
- /* for (j=1; j<=(cptcovs); j++) { /\* From model V1 + V2*age+ V3 + V3*V4 keeps V1 + V3 = 2 only *\/ */
- for (k=1; k<=cptcovt; k++) { /* From model V1 + V2*age + V3 + V3*V4 keeps V1 + V3 = 2 only */
- 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]) {
- case 0: /* Testing on fixed dummy covariate, simple or product of fixed */
- 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*/
- ij=(int)(covar[Tvar[k]][i]);
- /* ij=0 or 1 or -1. Value of the covariate Tvar[j] for individual i
- * If product of Vn*Vm, still boolean *:
- * If it was coded 1, 2, 3, 4 should be splitted into 3 boolean variables
- * 1 => 0 0 0, 2 => 0 0 1, 3 => 0 1 1, 4=1 0 0 */
- /* Finds for covariate j, n=Tvar[j] of Vn . ij is the
- modality of the nth covariate of individual i. */
- if (ij > modmaxcovj)
- modmaxcovj=ij;
- else if (ij < modmincovj)
- modmincovj=ij;
- if ((ij < -1) && (ij > NCOVMAX)){
- printf( "Error: minimal is less than -1 or maximal is bigger than %d. Exiting. \n", NCOVMAX );
- exit(1);
- }else
- Ndum[ij]++; /*counts and stores the occurence of this modality 0, 1, -1*/
- /* If coded 1, 2, 3 , counts the number of 1 Ndum[1], number of 2, Ndum[2], etc */
- /*printf("i=%d ij=%d Ndum[ij]=%d imx=%d",i,ij,Ndum[ij],imx);*/
- /* getting the maximum value of the modality of the covariate
- (should be 0 or 1 now) Tvar[j]. If V=sex and male is coded 0 and
- female ies 1, then modmaxcovj=1.
- */
- } /* end for loop on individuals i */
- printf(" Minimal and maximal values of %d th covariate V%d: min=%d max=%d \n", k, Tvar[k], modmincovj, modmaxcovj);
- fprintf(ficlog," Minimal and maximal values of %d th covariate V%d: min=%d max=%d \n", k, Tvar[k], modmincovj, modmaxcovj);
- cptcode=modmaxcovj;
- /* Ndum[0] = frequency of 0 for model-covariate j, Ndum[1] frequency of 1 etc. */
- /*for (i=0; i<=cptcode; i++) {*/
- for (j=modmincovj; j<=modmaxcovj; j++) { /* j=-1 ? 0 and 1*//* For each value j of the modality of model-cov k */
- printf("Frequencies of covariates %d ie V%d with value %d: %d\n", k, Tvar[k], j, Ndum[j]);
- fprintf(ficlog, "Frequencies of covariates %d ie V%d with value %d: %d\n", k, Tvar[k], j, Ndum[j]);
- if( Ndum[j] != 0 ){ /* Counts if nobody answered modality j ie empty modality, we skip it and reorder */
- if( j != -1){
- ncodemax[k]++; /* ncodemax[k]= Number of modalities of the k th
- covariate for which somebody answered excluding
- undefined. Usually 2: 0 and 1. */
- }
- ncodemaxwundef[k]++; /* ncodemax[j]= Number of modalities of the k th
- covariate for which somebody answered including
- undefined. Usually 3: -1, 0 and 1. */
- } /* In fact ncodemax[k]=2 (dichotom. variables only) but it could be more for
- * historical reasons: 3 if coded 1, 2, 3 and 4 and Ndum[2]=0 */
- } /* Ndum[-1] number of undefined modalities */
+ for (k=1; k <= maxncov; k++) ncodemax[k]=0; /* Horrible constant again replaced by NCOVMAX */
+
+ /* Loop on covariates without age and products and no quantitative variable */
+ /* for (j=1; j<=(cptcovs); j++) { /\* From model V1 + V2*age+ V3 + V3*V4 keeps V1 + V3 = 2 only *\/ */
+ for (k=1; k<=cptcovt; k++) { /* From model V1 + V2*age + V3 + V3*V4 keeps V1 + V3 = 2 only */
+ 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]) {
+ case 0: /* Testing on fixed dummy covariate, simple or product of fixed */
+ 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*/
+ ij=(int)(covar[Tvar[k]][i]);
+ /* ij=0 or 1 or -1. Value of the covariate Tvar[j] for individual i
+ * If product of Vn*Vm, still boolean *:
+ * If it was coded 1, 2, 3, 4 should be splitted into 3 boolean variables
+ * 1 => 0 0 0, 2 => 0 0 1, 3 => 0 1 1, 4=1 0 0 */
+ /* Finds for covariate j, n=Tvar[j] of Vn . ij is the
+ modality of the nth covariate of individual i. */
+ if (ij > modmaxcovj)
+ modmaxcovj=ij;
+ else if (ij < modmincovj)
+ modmincovj=ij;
+ if ((ij < -1) && (ij > NCOVMAX)){
+ printf( "Error: minimal is less than -1 or maximal is bigger than %d. Exiting. \n", NCOVMAX );
+ exit(1);
+ }else
+ Ndum[ij]++; /*counts and stores the occurence of this modality 0, 1, -1*/
+ /* If coded 1, 2, 3 , counts the number of 1 Ndum[1], number of 2, Ndum[2], etc */
+ /*printf("i=%d ij=%d Ndum[ij]=%d imx=%d",i,ij,Ndum[ij],imx);*/
+ /* getting the maximum value of the modality of the covariate
+ (should be 0 or 1 now) Tvar[j]. If V=sex and male is coded 0 and
+ female ies 1, then modmaxcovj=1.
+ */
+ } /* end for loop on individuals i */
+ printf(" Minimal and maximal values of %d th (fixed) covariate V%d: min=%d max=%d \n", k, Tvar[k], modmincovj, modmaxcovj);
+ fprintf(ficlog," Minimal and maximal values of %d th (fixed) covariate V%d: min=%d max=%d \n", k, Tvar[k], modmincovj, modmaxcovj);
+ cptcode=modmaxcovj;
+ /* Ndum[0] = frequency of 0 for model-covariate j, Ndum[1] frequency of 1 etc. */
+ /*for (i=0; i<=cptcode; i++) {*/
+ for (j=modmincovj; j<=modmaxcovj; j++) { /* j=-1 ? 0 and 1*//* For each value j of the modality of model-cov k */
+ printf("Frequencies of (fixed) covariate %d ie V%d with value %d: %d\n", k, Tvar[k], j, Ndum[j]);
+ fprintf(ficlog, "Frequencies of (fixed) covariate %d ie V%d with value %d: %d\n", k, Tvar[k], j, Ndum[j]);
+ if( Ndum[j] != 0 ){ /* Counts if nobody answered modality j ie empty modality, we skip it and reorder */
+ if( j != -1){
+ ncodemax[k]++; /* ncodemax[k]= Number of modalities of the k th
+ covariate for which somebody answered excluding
+ undefined. Usually 2: 0 and 1. */
+ }
+ ncodemaxwundef[k]++; /* ncodemax[j]= Number of modalities of the k th
+ covariate for which somebody answered including
+ undefined. Usually 3: -1, 0 and 1. */
+ } /* In fact ncodemax[k]=2 (dichotom. variables only) but it could be more for
+ * historical reasons: 3 if coded 1, 2, 3 and 4 and Ndum[2]=0 */
+ } /* Ndum[-1] number of undefined modalities */
- /* j is a covariate, n=Tvar[j] of Vn; Fills nbcode */
- /* For covariate j, modalities could be 1, 2, 3, 4, 5, 6, 7. */
- /* If Ndum[1]=0, Ndum[2]=0, Ndum[3]= 635, Ndum[4]=0, Ndum[5]=0, Ndum[6]=27, Ndum[7]=125; */
- /* modmincovj=3; modmaxcovj = 7; */
- /* There are only 3 modalities non empty 3, 6, 7 (or 2 if 27 is too few) : ncodemax[j]=3; */
- /* which will be coded 0, 1, 2 which in binary on 2=3-1 digits are 0=00 1=01, 2=10; */
- /* defining two dummy variables: variables V1_1 and V1_2.*/
- /* nbcode[Tvar[j]][ij]=k; */
- /* nbcode[Tvar[j]][1]=0; */
- /* nbcode[Tvar[j]][2]=1; */
- /* nbcode[Tvar[j]][3]=2; */
- /* To be continued (not working yet). */
- ij=0; /* ij is similar to i but can jump over null modalities */
- for (i=modmincovj; i<=modmaxcovj; i++) { /* i= 1 to 2 for dichotomous, or from 1 to 3 or from -1 or 0 to 1 currently*/
- if (Ndum[i] == 0) { /* If nobody responded to this modality k */
- break;
- }
- ij++;
- nbcode[Tvar[k]][ij]=i; /* stores the original value of modality i in an array nbcode, ij modality from 1 to last non-nul modality. nbcode[1][1]=0 nbcode[1][2]=1*/
- cptcode = ij; /* New max modality for covar j */
- } /* end of loop on modality i=-1 to 1 or more */
- break;
- case 1: /* Testing on varying covariate, could be simple and
- * should look at waves or product of fixed *
- * varying. No time to test -1, assuming 0 and 1 only */
- ij=0;
- for(i=0; i<=1;i++){
- nbcode[Tvar[k]][++ij]=i;
- }
- break;
- default:
- break;
- } /* end switch */
- } /* end dummy test */
+ /* j is a covariate, n=Tvar[j] of Vn; Fills nbcode */
+ /* For covariate j, modalities could be 1, 2, 3, 4, 5, 6, 7. */
+ /* If Ndum[1]=0, Ndum[2]=0, Ndum[3]= 635, Ndum[4]=0, Ndum[5]=0, Ndum[6]=27, Ndum[7]=125; */
+ /* modmincovj=3; modmaxcovj = 7; */
+ /* There are only 3 modalities non empty 3, 6, 7 (or 2 if 27 is too few) : ncodemax[j]=3; */
+ /* which will be coded 0, 1, 2 which in binary on 2=3-1 digits are 0=00 1=01, 2=10; */
+ /* defining two dummy variables: variables V1_1 and V1_2.*/
+ /* nbcode[Tvar[j]][ij]=k; */
+ /* nbcode[Tvar[j]][1]=0; */
+ /* nbcode[Tvar[j]][2]=1; */
+ /* nbcode[Tvar[j]][3]=2; */
+ /* To be continued (not working yet). */
+ ij=0; /* ij is similar to i but can jump over null modalities */
+ for (i=modmincovj; i<=modmaxcovj; i++) { /* i= 1 to 2 for dichotomous, or from 1 to 3 or from -1 or 0 to 1 currently*/
+ if (Ndum[i] == 0) { /* If nobody responded to this modality k */
+ break;
+ }
+ ij++;
+ nbcode[Tvar[k]][ij]=i; /* stores the original value of modality i in an array nbcode, ij modality from 1 to last non-nul modality. nbcode[1][1]=0 nbcode[1][2]=1*/
+ cptcode = ij; /* New max modality for covar j */
+ } /* end of loop on modality i=-1 to 1 or more */
+ break;
+ case 1: /* Testing on varying covariate, could be simple and
+ * should look at waves or product of fixed *
+ * varying. No time to test -1, assuming 0 and 1 only */
+ ij=0;
+ for(i=0; i<=1;i++){
+ nbcode[Tvar[k]][++ij]=i;
+ }
+ break;
+ default:
+ break;
+ } /* end switch */
+ } /* end dummy test */
- /* for (k=0; k<= cptcode; k++) { /\* k=-1 ? k=0 to 1 *\//\* Could be 1 to 4 *\//\* cptcode=modmaxcovj *\/ */
- /* /\*recode from 0 *\/ */
- /* k is a modality. If we have model=V1+V1*sex */
- /* then: nbcode[1][1]=0 ; nbcode[1][2]=1; nbcode[2][1]=0 ; nbcode[2][2]=1; */
- /* But if some modality were not used, it is recoded from 0 to a newer modmaxcovj=cptcode *\/ */
- /* } */
- /* /\* cptcode = ij; *\/ /\* New max modality for covar j *\/ */
- /* if (ij > ncodemax[j]) { */
- /* printf( " Error ij=%d > ncodemax[%d]=%d\n", ij, j, ncodemax[j]); */
- /* fprintf(ficlog, " Error ij=%d > ncodemax[%d]=%d\n", ij, j, ncodemax[j]); */
- /* break; */
- /* } */
- /* } /\* end of loop on modality k *\/ */
- } /* end of loop on model-covariate j. nbcode[Tvarj][1]=0 and nbcode[Tvarj][2]=1 sets the value of covariate j*/
-
- for (k=-1; k< maxncov; k++) Ndum[k]=0;
- /* Look at fixed dummy (single or product) covariates to check empty modalities */
- for (i=1; i<=ncovmodel-2-nagesqr; i++) { /* -2, cste and age and eventually age*age */
- /* Listing of all covariables in statement model to see if some covariates appear twice. For example, V1 appears twice in V1+V1*V2.*/
- ij=Tvar[i]; /* Tvar 5,4,3,6,5,7,1,4 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V4*age */
- Ndum[ij]++; /* Count the # of 1, 2 etc: {1,1,1,2,2,1,1} because V1 once, V2 once, two V4 and V5 in above */
- /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, {2, 1, 1, 1, 2, 1, 1, 0, 0} */
- } /* V4+V3+V5, Ndum[1]@5={0, 0, 1, 1, 1} */
-
- 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) */
- /*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 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, */
- 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 */
- if(Fixed[k]!=0)
- anyvaryingduminmodel=1;
- /* }else if((Ndum[i]!=0) && (i<=ncovcol+nqv)){ */
- /* Tvaraff[++ij]=-10; /\* Dont'n know how to treat quantitative variables yet *\/ */
- /* }else if((Ndum[i]!=0) && (i<=ncovcol+nqv+ntv)){ */
- /* Tvaraff[++ij]=i; /\*For printing (unclear) *\/ */
- /* }else if((Ndum[i]!=0) && (i<=ncovcol+nqv+ntv+nqtv)){ */
- /* Tvaraff[++ij]=-20; /\* Dont'n know how to treat quantitative variables yet *\/ */
- }
- } /* Tvaraff[1]@5 {3, 4, -20, 0, 0} Very strange */
- /* ij--; */
- /* cptcoveff=ij; /\*Number of total covariates*\/ */
- *cptcov=ij; /*Number of total real effective covariates: effective
- * 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++){
- Tvaraff[j]=0;
- Tmodelind[j]=0;
- }
- for(j=ntveff+1; j<= cptcovt; j++){
- TmodelInvind[j]=0;
- }
- /* To be sorted */
- ;
-}
+ /* for (k=0; k<= cptcode; k++) { /\* k=-1 ? k=0 to 1 *\//\* Could be 1 to 4 *\//\* cptcode=modmaxcovj *\/ */
+ /* /\*recode from 0 *\/ */
+ /* k is a modality. If we have model=V1+V1*sex */
+ /* then: nbcode[1][1]=0 ; nbcode[1][2]=1; nbcode[2][1]=0 ; nbcode[2][2]=1; */
+ /* But if some modality were not used, it is recoded from 0 to a newer modmaxcovj=cptcode *\/ */
+ /* } */
+ /* /\* cptcode = ij; *\/ /\* New max modality for covar j *\/ */
+ /* if (ij > ncodemax[j]) { */
+ /* printf( " Error ij=%d > ncodemax[%d]=%d\n", ij, j, ncodemax[j]); */
+ /* fprintf(ficlog, " Error ij=%d > ncodemax[%d]=%d\n", ij, j, ncodemax[j]); */
+ /* break; */
+ /* } */
+ /* } /\* end of loop on modality k *\/ */
+ } /* end of loop on model-covariate j. nbcode[Tvarj][1]=0 and nbcode[Tvarj][2]=1 sets the value of covariate j*/
+
+ for (k=-1; k< maxncov; k++) Ndum[k]=0;
+ /* Look at fixed dummy (single or product) covariates to check empty modalities */
+ for (i=1; i<=ncovmodel-2-nagesqr; i++) { /* -2, cste and age and eventually age*age */
+ /* Listing of all covariables in statement model to see if some covariates appear twice. For example, V1 appears twice in V1+V1*V2.*/
+ ij=Tvar[i]; /* Tvar 5,4,3,6,5,7,1,4 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V4*age */
+ Ndum[ij]++; /* Count the # of 1, 2 etc: {1,1,1,2,2,1,1} because V1 once, V2 once, two V4 and V5 in above */
+ /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, {2, 1, 1, 1, 2, 1, 1, 0, 0} */
+ } /* V4+V3+V5, Ndum[1]@5={0, 0, 1, 1, 1} */
+
+ 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) */
+ /*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 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, */
+ 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 */
+ if(Fixed[k]!=0)
+ anyvaryingduminmodel=1;
+ /* }else if((Ndum[i]!=0) && (i<=ncovcol+nqv)){ */
+ /* Tvaraff[++ij]=-10; /\* Dont'n know how to treat quantitative variables yet *\/ */
+ /* }else if((Ndum[i]!=0) && (i<=ncovcol+nqv+ntv)){ */
+ /* Tvaraff[++ij]=i; /\*For printing (unclear) *\/ */
+ /* }else if((Ndum[i]!=0) && (i<=ncovcol+nqv+ntv+nqtv)){ */
+ /* Tvaraff[++ij]=-20; /\* Dont'n know how to treat quantitative variables yet *\/ */
+ }
+ } /* Tvaraff[1]@5 {3, 4, -20, 0, 0} Very strange */
+ /* ij--; */
+ /* cptcoveff=ij; /\*Number of total covariates*\/ */
+ *cptcov=ij; /*Number of total real effective covariates: effective
+ * 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++){
+ Tvaraff[j]=0;
+ Tmodelind[j]=0;
+ }
+ for(j=ntveff+1; j<= cptcovt; j++){
+ TmodelInvind[j]=0;
+ }
+ /* To be sorted */
+ ;
+ }
/*********** Health Expectancies ****************/
xp[i] = x[i] + (i==theta ?delti[theta]:0);
}
- prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij, nresult);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij, nres);
if (popbased==1) {
if(mobilav ==0){
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, nresult);
+ prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp, ij, nres);
if (popbased==1) {
if(mobilav ==0){
/* end ppptj */
/* x centered again */
- prevalim(prlim,nlstate,x,age,oldm,savm,ftolpl,ncvyearp,ij, nresult);
+ prevalim(prlim,nlstate,x,age,oldm,savm,ftolpl,ncvyearp,ij, nres);
if (popbased==1) {
if(mobilav ==0){
if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
else fprintf(ficgp," %%*lf (%%*lf)");
}
- fprintf(ficgp,"\" t\"Period (stable) prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2==%d ? $4+1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1,nres);
+ fprintf(ficgp,"\" t\"Period (stable) prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2==%d ? $3+1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1,nres);
for (i=1; i<= nlstate ; i ++) {
if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
else fprintf(ficgp," %%*lf (%%*lf)");
}
- fprintf(ficgp,"\" t\"95%% CI\" w l lt 1,\"%s\" every :::%d::%d u 1:($2==%d ? $4-1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1,nres);
+ fprintf(ficgp,"\" t\"95%% CI\" w l lt 1,\"%s\" every :::%d::%d u 1:($2==%d ? $3-1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1,nres);
for (i=1; i<= nlstate ; i ++) {
if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
else fprintf(ficgp," %%*lf (%%*lf)");
fprintf(ficgp,"\" t\"\" w l lt 1,\"%s\" every :::%d::%d u 1:($%d) t\"Observed prevalence\" w l lt 2",subdirf2(fileresu,"P_"),k1-1,k1-1,2+4*(cpt-1));
if(backcast==1){ /* We need to get the corresponding values of the covariates involved in this combination k1 */
/* fprintf(ficgp,",\"%s\" every :::%d::%d u 1:($%d) t\"Backward stable prevalence\" w l lt 3",subdirf2(fileresu,"PLB_"),k1-1,k1-1,1+cpt); */
- fprintf(ficgp,",\"%s\" u ($2==%d ?$1:1/0):(",subdirf2(fileresu,"PLB_"),nres); /* Age is in 1, nres in 2 to be fixed */
+ fprintf(ficgp,",\"%s\" u 1:((",subdirf2(fileresu,"PLB_")); /* Age is in 1, nres in 2 to be fixed */
if(cptcoveff ==0){
fprintf(ficgp,"$%d)) t 'Backward prevalence in state %d' with line ", 2+(cpt-1), cpt );
}else{
/* '' u 6:(($1==1 && $2==0 && $3==2 && $4==0)? $9/(1.-$15) : 1/0):($5==2000? 3:2) t 'p.1' with line lc variable*/
if(k==cptcoveff){
fprintf(ficgp,"$%d==%d && $%d==%d)? $%d : 1/0) t 'Backward prevalence in state %d' ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv], \
- 4+(cpt-1), cpt ); /* 4 or 6 ?*/
+ 2+cptcoveff*2+(cpt-1), cpt ); /* 4 or 6 ?*/
}else{
fprintf(ficgp,"$%d==%d && $%d==%d && ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv]);
kl++;
TvarFind[ncovf]=k;
TvarFQ[nqfveff]=Tvar[k]-ncovcol; /* TvarFQ[1]=V2-1=1st in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
TvarFQind[nqfveff]=k; /* TvarFQind[1]=6 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
- }else if( Tvar[k] <=ncovcol+nqv+ntv && Typevar[k]==0){/* Only simple time varying variables */
+ }else if( Tvar[k] <=ncovcol+nqv+ntv && Typevar[k]==0){/* Only simple time varying dummy variables */
Fixed[k]= 1;
Dummy[k]= 0;
ntveff++; /* Only simple time varying dummy variable */
TvarsDind[nsd]=k;
ncovv++; /* Only simple time varying variables */
TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=k;
+ TvarVind[ncovv]=k; /* TvarVind[2]=2 TvarVind[3]=3 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Any time varying singele */
TvarVD[ntveff]=Tvar[k]; /* TvarVD[1]=V4 TvarVD[2]=V3 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple time varying dummy variable */
TvarVDind[ntveff]=k; /* TvarVDind[1]=2 TvarVDind[2]=3 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple time varying dummy variable */
printf("Quasi Tmodelind[%d]=%d,Tvar[Tmodelind[%d]]=V%d, ncovcol=%d, nqv=%d,Tvar[k]- ncovcol-nqv=%d\n",ntveff,k,ntveff,Tvar[k], ncovcol, nqv,Tvar[k]- ncovcol-nqv);
TvarsQ[nsq]=Tvar[k];
TvarsQind[nsq]=k;
TvarV[ncovv]=Tvar[k];
- TvarVind[ncovv]=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 */
TvarVQind[nqtveff]=k; /* TvarVQind[1]=1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple time varying quantitative variable */
TmodelInvQind[nqtveff]=Tvar[k]- ncovcol-nqv-ntv;/* Only simple time varying quantitative variable */
if(mobilavproj > 0){
/* bprevalim(bprlim, mobaverage, nlstate, p, age, ageminpar, agemaxpar, oldm, savm, doldm, dsavm, ftolpl, ncvyearp, k); */
/* bprevalim(bprlim, mobaverage, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
- bprevalim(bprlim, mobaverage, nlstate, p, age, ftolpl, ncvyearp, k);
+ bprevalim(bprlim, mobaverage, nlstate, p, age, ftolpl, ncvyearp, k, nres);
}else if (mobilavproj == 0){
printf("There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
fprintf(ficlog,"There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
exit(1);
}else{
/* bprevalim(bprlim, probs, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
- bprevalim(bprlim, probs, nlstate, p, age, ftolpl, ncvyearp, k);
+ bprevalim(bprlim, probs, nlstate, p, age, ftolpl, ncvyearp, k,nres);
}
fprintf(ficresplb,"%.0f ",age );
for(j=1;j<=cptcoveff;j++)
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:");
+ printf("\n#****** Result for:");
+ fprintf(ficrest,"\n#****** Result for:");
+ fprintf(ficlog,"\n#****** Result for:");
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)]);
git://henry.ined.fr / .git / commitdiff