/* $Id$
$State$
$Log$
+ Revision 1.230 2016/08/22 06:55:53 brouard
+ Summary: Not working
+
Revision 1.229 2016/07/23 09:45:53 brouard
Summary: Completing for func too
double ***cotqvar; /* Time varying quantitative covariate itqv */
double idx;
int **nbcode, *Tvar; /**< model=V2 => Tvar[1]= 2 */
+int *TvarFD; /**< TvarFD[1]=V1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
+int *TvarFDind; /* TvarFDind[1]=9 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
+int *TvarFQ; /* TvarFQ[1]=V2 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
+int *TvarFQind; /* TvarFQind[1]=6 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
+int *TvarVD; /* TvarVD[1]=V5 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
+int *TvarVDind; /* TvarVDind[1]=1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
+int *TvarVQ; /* TvarVQ[1]=V5 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple time varying quantitative variable */
+int *TvarVQind; /* TvarVQind[1]=1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple time varying quantitative variable */
+
int *Tvarsel; /**< Selected covariates for output */
double *Tvalsel; /**< Selected modality value of covariate for output */
int *Typevar; /**< 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product */
int cptcovprod, *Tvaraff, *invalidvarcomb;
double *lsurv, *lpop, *tpop;
+#define FD 1; /* Fixed dummy covariate */
+#define FQ 2; /* Fixed quantitative covariate */
+#define FP 3; /* Fixed product covariate */
+#define FPDD 7; /* Fixed product dummy*dummy covariate */
+#define FPDQ 8; /* Fixed product dummy*quantitative covariate */
+#define FPQQ 9; /* Fixed product quantitative*quantitative covariate */
+#define VD 10; /* Varying dummy covariate */
+#define VQ 11; /* Varying quantitative covariate */
+#define VP 12; /* Varying product covariate */
+#define VPDD 13; /* Varying product dummy*dummy covariate */
+#define VPDQ 14; /* Varying product dummy*quantitative covariate */
+#define VPQQ 15; /* Varying product quantitative*quantitative covariate */
+#define APFD 16; /* Age product * fixed dummy covariate */
+#define APFQ 17; /* Age product * fixed quantitative covariate */
+#define APVD 18; /* Age product * varying dummy covariate */
+#define APVQ 19; /* Age product * varying quantitative covariate */
+
+#define FTYPE 1; /* Fixed covariate */
+#define VTYPE 2; /* Varying covariate (loop in wave) */
+#define ATYPE 2; /* Age product covariate (loop in dh within wave)*/
+
+struct kmodel{
+ int maintype; /* main type */
+ int subtype; /* subtype */
+};
+struct kmodel modell[NCOVMAX];
+
double ftol=FTOL; /**< Tolerance for computing Max Likelihood */
double ftolhess; /**< Tolerance for computing hessian */
meaning that decodemodel should be used cotvar[mw[mi+1][i]][TTvar[iv]][i]
*/
for(mi=1; mi<= wav[i]-1; mi++){
- for(itv=1; itv <= ntveff; itv++){ /* Varying dummy covariates */
- /* cov[ioffset+itv]=cotvar[mw[mi][i]][Tvar[itv]][i]; /\* Not sure, Tvar V4+V3+V5 Tvaraff ? *\/ */
- cov[ioffset+itv]=cotvar[mw[mi][i]][TmodelInvind[itv]][i];
- }
- for(iqtv=1; iqtv <= nqtveff; iqtv++){ /* Varying quantitatives covariates */
- if(cotqvar[mw[mi][i]][iqtv][i] == -1){
- printf("i=%d, mi=%d, iqtv=%d, cotqvar[mw[mi][i]][iqtv][i]=%f",i,mi,iqtv,cotqvar[mw[mi][i]][iqtv][i]);
- }
- cov[ioffset+ntveff+iqtv]=cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i];
- /* cov[ioffset+ntveff+iqtv]=cotqvar[mw[mi][i]][iqtv][i]; */
- }
- /* ioffset=2+nagesqr+cptcovn+nqv+ntv+nqtv; */
- 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(d=0; d<dh[mi][i]; d++){
- newm=savm;
- agexact=agev[mw[mi][i]][i]+d*stepm/YEARM;
- cov[2]=agexact;
- if(nagesqr==1)
- cov[3]= agexact*agexact; /* Should be changed here */
- for (kk=1; kk<=cptcovage;kk++) {
- cov[Tage[kk]+2+nagesqr]=covar[Tvar[Tage[kk]]][i]*agexact; /* Tage[kk] gives the data-covariate associated with age */
- }
- 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 */
+ for(itv=1; itv <= ntveff; itv++){ /* Varying dummy covariates */
+ /* cov[ioffset+itv]=cotvar[mw[mi][i]][Tvar[itv]][i]; /\* Not sure, Tvar V4+V3+V5 Tvaraff ? *\/ */
+ cov[ioffset+itv]=cotvar[mw[mi][i]][TmodelInvind[itv]][i];
+ }
+ for(iqtv=1; iqtv <= nqtveff; iqtv++){ /* Varying quantitatives covariates */
+ if(cotqvar[mw[mi][i]][iqtv][i] == -1){
+ printf("i=%d, mi=%d, iqtv=%d, cotqvar[mw[mi][i]][iqtv][i]=%f",i,mi,iqtv,cotqvar[mw[mi][i]][iqtv][i]);
+ }
+ cov[ioffset+ntveff+iqtv]=cotqvar[mw[mi][i]][TmodelInvQind[iqtv]][i];
+ /* cov[ioffset+ntveff+iqtv]=cotqvar[mw[mi][i]][iqtv][i]; */
+ }
+ /* ioffset=2+nagesqr+cptcovn+nqv+ntv+nqtv; */
+ 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(d=0; d<dh[mi][i]; d++){
+ newm=savm;
+ agexact=agev[mw[mi][i]][i]+d*stepm/YEARM;
+ cov[2]=agexact;
+ if(nagesqr==1)
+ cov[3]= agexact*agexact; /* Should be changed here */
+ for (kk=1; kk<=cptcovage;kk++) {
+ cov[Tage[kk]+2+nagesqr]=covar[Tvar[Tage[kk]]][i]*agexact; /* Tage[kk] gives the data-covariate associated with age */
+ }
+ 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 */
- /*lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/ /* Original formula */
- /* But now since version 0.9 we anticipate for bias at large stepm.
- * If stepm is larger than one month (smallest stepm) and if the exact delay
- * (in months) between two waves is not a multiple of stepm, we rounded to
- * the nearest (and in case of equal distance, to the lowest) interval but now
- * we keep into memory the bias bh[mi][i] and also the previous matrix product
- * (i.e to dh[mi][i]-1) saved in 'savm'. Then we inter(extra)polate the
- * probability in order to take into account the bias as a fraction of the way
- * from savm to out if bh is negative or even beyond if bh is positive. bh varies
- * -stepm/2 to stepm/2 .
- * For stepm=1 the results are the same as for previous versions of Imach.
- * For stepm > 1 the results are less biased than in previous versions.
- */
- s1=s[mw[mi][i]][i];
- s2=s[mw[mi+1][i]][i];
- bbh=(double)bh[mi][i]/(double)stepm;
- /* bias bh is positive if real duration
- * is higher than the multiple of stepm and negative otherwise.
- */
- /* lli= (savm[s1][s2]>1.e-8 ?(1.+bbh)*log(out[s1][s2])- bbh*log(savm[s1][s2]):log((1.+bbh)*out[s1][s2]));*/
- if( s2 > nlstate){
- /* i.e. if s2 is a death state and if the date of death is known
- then the contribution to the likelihood is the probability to
- die between last step unit time and current step unit time,
- which is also equal to probability to die before dh
- minus probability to die before dh-stepm .
- In version up to 0.92 likelihood was computed
- as if date of death was unknown. Death was treated as any other
- health state: the date of the interview describes the actual state
- and not the date of a change in health state. The former idea was
- to consider that at each interview the state was recorded
- (healthy, disable or death) and IMaCh was corrected; but when we
- introduced the exact date of death then we should have modified
- the contribution of an exact death to the likelihood. This new
- contribution is smaller and very dependent of the step unit
- stepm. It is no more the probability to die between last interview
- and month of death but the probability to survive from last
- interview up to one month before death multiplied by the
- probability to die within a month. Thanks to Chris
- Jackson for correcting this bug. Former versions increased
- mortality artificially. The bad side is that we add another loop
- which slows down the processing. The difference can be up to 10%
- lower mortality.
- */
- /* If, at the beginning of the maximization mostly, the
- cumulative probability or probability to be dead is
- constant (ie = 1) over time d, the difference is equal to
- 0. out[s1][3] = savm[s1][3]: probability, being at state
- s1 at precedent wave, to be dead a month before current
- wave is equal to probability, being at state s1 at
- precedent wave, to be dead at mont of the current
- wave. Then the observed probability (that this person died)
- is null according to current estimated parameter. In fact,
- it should be very low but not zero otherwise the log go to
- infinity.
- */
+ /*lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/ /* Original formula */
+ /* But now since version 0.9 we anticipate for bias at large stepm.
+ * If stepm is larger than one month (smallest stepm) and if the exact delay
+ * (in months) between two waves is not a multiple of stepm, we rounded to
+ * the nearest (and in case of equal distance, to the lowest) interval but now
+ * we keep into memory the bias bh[mi][i] and also the previous matrix product
+ * (i.e to dh[mi][i]-1) saved in 'savm'. Then we inter(extra)polate the
+ * probability in order to take into account the bias as a fraction of the way
+ * from savm to out if bh is negative or even beyond if bh is positive. bh varies
+ * -stepm/2 to stepm/2 .
+ * For stepm=1 the results are the same as for previous versions of Imach.
+ * For stepm > 1 the results are less biased than in previous versions.
+ */
+ s1=s[mw[mi][i]][i];
+ s2=s[mw[mi+1][i]][i];
+ bbh=(double)bh[mi][i]/(double)stepm;
+ /* bias bh is positive if real duration
+ * is higher than the multiple of stepm and negative otherwise.
+ */
+ /* lli= (savm[s1][s2]>1.e-8 ?(1.+bbh)*log(out[s1][s2])- bbh*log(savm[s1][s2]):log((1.+bbh)*out[s1][s2]));*/
+ if( s2 > nlstate){
+ /* i.e. if s2 is a death state and if the date of death is known
+ then the contribution to the likelihood is the probability to
+ die between last step unit time and current step unit time,
+ which is also equal to probability to die before dh
+ minus probability to die before dh-stepm .
+ In version up to 0.92 likelihood was computed
+ as if date of death was unknown. Death was treated as any other
+ health state: the date of the interview describes the actual state
+ and not the date of a change in health state. The former idea was
+ to consider that at each interview the state was recorded
+ (healthy, disable or death) and IMaCh was corrected; but when we
+ introduced the exact date of death then we should have modified
+ the contribution of an exact death to the likelihood. This new
+ contribution is smaller and very dependent of the step unit
+ stepm. It is no more the probability to die between last interview
+ and month of death but the probability to survive from last
+ interview up to one month before death multiplied by the
+ probability to die within a month. Thanks to Chris
+ Jackson for correcting this bug. Former versions increased
+ mortality artificially. The bad side is that we add another loop
+ which slows down the processing. The difference can be up to 10%
+ lower mortality.
+ */
+ /* If, at the beginning of the maximization mostly, the
+ cumulative probability or probability to be dead is
+ constant (ie = 1) over time d, the difference is equal to
+ 0. out[s1][3] = savm[s1][3]: probability, being at state
+ s1 at precedent wave, to be dead a month before current
+ wave is equal to probability, being at state s1 at
+ precedent wave, to be dead at mont of the current
+ wave. Then the observed probability (that this person died)
+ is null according to current estimated parameter. In fact,
+ it should be very low but not zero otherwise the log go to
+ infinity.
+ */
/* #ifdef INFINITYORIGINAL */
/* lli=log(out[s1][s2] - savm[s1][s2]); */
/* #else */
for (i=1,ipmx=0, sw=0.; i<=imx; i++){
ioffset=2+nagesqr+cptcovage;
/* for (k=1; k<=cptcovn;k++) cov[2+nagesqr+k]=covar[Tvar[k]][i]; */
- for (k=1; k<=ncoveff+nqfveff;k++){ /* Simple and product fixed Dummy covariates without age* products */
+ for (k=1; k<=ncoveff;k++){ /* Simple and product fixed Dummy covariates without age* products */
cov[++ioffset]=covar[TvarFD[k]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V1 is fixed (k=6)*/
}
- for(iqv=1; iqv <= nqfveff; iqv++){ /* Quantitative fixed covariates */
- cov[++ioffset]=coqvar[Tvar[iqv]][i]; /* Only V1 k=9 */
+ for (k=1; k<=nqfveff;k++){ /* Simple and product fixed Quantitative covariates without age* products */
+ cov[++ioffset]=coqvar[TvarFQ[k]][i];/* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, only V2 and V1*V2 is fixed (k=6 and 7?)*/
}
+ /* for(iqv=1; iqv <= nqfveff; iqv++){ /\* Quantitative fixed covariates *\/ */
+ /* cov[++ioffset]=coqvar[Tvar[iqv]][i]; /\* Only V2 k=6 and V1*V2 7 *\/ */
+ /* } */
for(mi=1; mi<= wav[i]-1; mi++){ /* Varying with waves */
- for(itv=1; itv <= ntveff; itv++){ /* Varying dummy covariates */
- /* 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(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]); */
}
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++) {
+ 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;
} /* end mult */
s1=s[mw[mi][i]][i];
scanf("%d", i);*/
for (i=-5; i<=nlstate+ndeath; i++)
for (jk=-5; jk<=nlstate+ndeath; jk++)
- for(m=iagemin; m <= iagemax+3; m++)
- freq[i][jk][m]=0;
-
+ for(m=iagemin; m <= iagemax+3; m++)
+ freq[i][jk][m]=0;
+
for (i=1; i<=nlstate; i++) {
for(m=iagemin; m <= iagemax+3; m++)
- prop[i][m]=0;
+ prop[i][m]=0;
posprop[i]=0;
pospropt[i]=0;
}
/* meanqt[m][z1]=0.; */
/* } */
/* } */
-
+
dateintsum=0;
k2cpt=0;
/* For that combination of covariate j1, we count and print the frequencies in one pass */
for (iind=1; iind<=imx; iind++) { /* For each individual iind */
bool=1;
if(anyvaryingduminmodel==0){ /* If All fixed covariates */
- if (cptcoveff >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
+ if (cptcoveff >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
/* for (z1=1; z1<= nqfveff; z1++) { */
/* meanq[z1]+=coqvar[Tvar[z1]][iind]; /\* Computes mean of quantitative with selected filter *\/ */
/* } */
- for (z1=1; z1<=cptcoveff; z1++) {
- /* if(Tvaraff[z1] ==-20){ */
- /* /\* sumnew+=cotvar[mw[mi][iind]][z1][iind]; *\/ */
- /* }else if(Tvaraff[z1] ==-10){ */
- /* /\* sumnew+=coqvar[z1][iind]; *\/ */
- /* }else */
- if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]){
- /* Tests if this individual iind responded to j1 (V4=1 V3=0) */
- bool=0;
- /* printf("bool=%d i=%d, z1=%d, Tvaraff[%d]=%d, covar[Tvarff][%d]=%2f, codtabm(%d,%d)=%d, nbcode[Tvaraff][codtabm(%d,%d)=%d, j1=%d\n",
- bool,i,z1, z1, Tvaraff[z1],i,covar[Tvaraff[z1]][i],j1,z1,codtabm(j1,z1),
- j1,z1,nbcode[Tvaraff[z1]][codtabm(j1,z1)],j1);*/
- /* For j1=7 in V1+V2+V3+V4 = 0 1 1 0 and codtabm(7,3)=1 and nbcde[3][?]=1*/
- } /* Onlyf fixed */
- } /* end z1 */
- } /* cptcovn > 0 */
+ for (z1=1; z1<=cptcoveff; z1++) {
+ /* if(Tvaraff[z1] ==-20){ */
+ /* /\* sumnew+=cotvar[mw[mi][iind]][z1][iind]; *\/ */
+ /* }else if(Tvaraff[z1] ==-10){ */
+ /* /\* sumnew+=coqvar[z1][iind]; *\/ */
+ /* }else */
+ if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]){
+ /* Tests if this individual iind responded to j1 (V4=1 V3=0) */
+ bool=0;
+ /* printf("bool=%d i=%d, z1=%d, Tvaraff[%d]=%d, covar[Tvarff][%d]=%2f, codtabm(%d,%d)=%d, nbcode[Tvaraff][codtabm(%d,%d)=%d, j1=%d\n",
+ bool,i,z1, z1, Tvaraff[z1],i,covar[Tvaraff[z1]][i],j1,z1,codtabm(j1,z1),
+ j1,z1,nbcode[Tvaraff[z1]][codtabm(j1,z1)],j1);*/
+ /* For j1=7 in V1+V2+V3+V4 = 0 1 1 0 and codtabm(7,3)=1 and nbcde[3][?]=1*/
+ } /* Onlyf fixed */
+ } /* end z1 */
+ } /* cptcovn > 0 */
} /* end any */
if (bool==1){ /* We selected an individual iind satisfying combination j1 or all fixed */
- /* for(m=firstpass; m<=lastpass; m++){ */
- for(mi=1; mi<wav[iind];mi++){ /* For that wave */
- m=mw[mi][iind];
- if(anyvaryingduminmodel==1){ /* Some are varying covariates */
- for (z1=1; z1<=cptcoveff; z1++) {
- if( Fixed[Tmodelind[z1]]==1){
- iv= Tvar[Tmodelind[z1]]-ncovcol-nqv;
- if (cotvar[m][iv][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]) /* iv=1 to ntv, right modality */
- bool=0;
- }else if( Fixed[Tmodelind[z1]]== 0) { /* fixed */
- if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]) {
- bool=0;
- }
- }
- }
- }/* Some are varying covariates, we tried to speed up if all fixed covariates in the model, avoiding waves loop */
- /* bool =0 we keep that guy which corresponds to the combination of dummy values */
- if(bool==1){
- /* dh[m][iind] or dh[mw[mi][iind]][iind] is the delay between two effective (mi) waves m=mw[mi][iind]
- and mw[mi+1][iind]. dh depends on stepm. */
- agebegin=agev[m][iind]; /* Age at beginning of wave before transition*/
- ageend=agev[m][iind]+(dh[m][iind])*stepm/YEARM; /* Age at end of wave and transition */
- if(m >=firstpass && m <=lastpass){
- k2=anint[m][iind]+(mint[m][iind]/12.);
- /*if ((k2>=dateprev1) && (k2<=dateprev2)) {*/
- if(agev[m][iind]==0) agev[m][iind]=iagemax+1; /* All ages equal to 0 are in iagemax+1 */
- if(agev[m][iind]==1) agev[m][iind]=iagemax+2; /* All ages equal to 1 are in iagemax+2 */
- if (s[m][iind]>0 && s[m][iind]<=nlstate) /* If status at wave m is known and a live state */
- prop[s[m][iind]][(int)agev[m][iind]] += weight[iind]; /* At age of beginning of transition, where status is known */
- if (m<lastpass) {
- /* if(s[m][iind]==4 && s[m+1][iind]==4) */
- /* printf(" num=%ld m=%d, iind=%d s1=%d s2=%d agev at m=%d\n", num[iind], m, iind,s[m][iind],s[m+1][iind], (int)agev[m][iind]); */
- if(s[m][iind]==-1)
- printf(" num=%ld m=%d, iind=%d s1=%d s2=%d agev at m=%d agebegin=%.2f ageend=%.2f, agemed=%d\n", num[iind], m, iind,s[m][iind],s[m+1][iind], (int)agev[m][iind],agebegin, ageend, (int)((agebegin+ageend)/2.));
- freq[s[m][iind]][s[m+1][iind]][(int)agev[m][iind]] += weight[iind]; /* At age of beginning of transition, where status is known */
- /* freq[s[m][iind]][s[m+1][iind]][(int)((agebegin+ageend)/2.)] += weight[iind]; */
- freq[s[m][iind]][s[m+1][iind]][iagemax+3] += weight[iind]; /* Total is in iagemax+3 *//* At age of beginning of transition, where status is known */
- }
- } /* end if between passes */
- if ((agev[m][iind]>1) && (agev[m][iind]< (iagemax+3)) && (anint[m][iind]!=9999) && (mint[m][iind]!=99)) {
- dateintsum=dateintsum+k2;
- k2cpt++;
- /* printf("iind=%ld dateintmean = %lf dateintsum=%lf k2cpt=%lf k2=%lf\n",iind, dateintsum/k2cpt, dateintsum,k2cpt, k2); */
- }
- } /* end bool 2 */
- } /* end m */
+ /* for(m=firstpass; m<=lastpass; m++){ */
+ for(mi=1; mi<wav[iind];mi++){ /* For that wave */
+ m=mw[mi][iind];
+ if(anyvaryingduminmodel==1){ /* Some are varying covariates */
+ for (z1=1; z1<=cptcoveff; z1++) {
+ if( Fixed[Tmodelind[z1]]==1){
+ iv= Tvar[Tmodelind[z1]]-ncovcol-nqv;
+ if (cotvar[m][iv][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]) /* iv=1 to ntv, right modality */
+ bool=0;
+ }else if( Fixed[Tmodelind[z1]]== 0) { /* fixed */
+ if (covar[Tvaraff[z1]][iind]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]) {
+ bool=0;
+ }
+ }
+ }
+ }/* Some are varying covariates, we tried to speed up if all fixed covariates in the model, avoiding waves loop */
+ /* bool =0 we keep that guy which corresponds to the combination of dummy values */
+ if(bool==1){
+ /* dh[m][iind] or dh[mw[mi][iind]][iind] is the delay between two effective (mi) waves m=mw[mi][iind]
+ and mw[mi+1][iind]. dh depends on stepm. */
+ agebegin=agev[m][iind]; /* Age at beginning of wave before transition*/
+ ageend=agev[m][iind]+(dh[m][iind])*stepm/YEARM; /* Age at end of wave and transition */
+ if(m >=firstpass && m <=lastpass){
+ k2=anint[m][iind]+(mint[m][iind]/12.);
+ /*if ((k2>=dateprev1) && (k2<=dateprev2)) {*/
+ if(agev[m][iind]==0) agev[m][iind]=iagemax+1; /* All ages equal to 0 are in iagemax+1 */
+ if(agev[m][iind]==1) agev[m][iind]=iagemax+2; /* All ages equal to 1 are in iagemax+2 */
+ if (s[m][iind]>0 && s[m][iind]<=nlstate) /* If status at wave m is known and a live state */
+ prop[s[m][iind]][(int)agev[m][iind]] += weight[iind]; /* At age of beginning of transition, where status is known */
+ if (m<lastpass) {
+ /* if(s[m][iind]==4 && s[m+1][iind]==4) */
+ /* printf(" num=%ld m=%d, iind=%d s1=%d s2=%d agev at m=%d\n", num[iind], m, iind,s[m][iind],s[m+1][iind], (int)agev[m][iind]); */
+ if(s[m][iind]==-1)
+ printf(" num=%ld m=%d, iind=%d s1=%d s2=%d agev at m=%d agebegin=%.2f ageend=%.2f, agemed=%d\n", num[iind], m, iind,s[m][iind],s[m+1][iind], (int)agev[m][iind],agebegin, ageend, (int)((agebegin+ageend)/2.));
+ freq[s[m][iind]][s[m+1][iind]][(int)agev[m][iind]] += weight[iind]; /* At age of beginning of transition, where status is known */
+ /* freq[s[m][iind]][s[m+1][iind]][(int)((agebegin+ageend)/2.)] += weight[iind]; */
+ freq[s[m][iind]][s[m+1][iind]][iagemax+3] += weight[iind]; /* Total is in iagemax+3 *//* At age of beginning of transition, where status is known */
+ }
+ } /* end if between passes */
+ if ((agev[m][iind]>1) && (agev[m][iind]< (iagemax+3)) && (anint[m][iind]!=9999) && (mint[m][iind]!=99)) {
+ dateintsum=dateintsum+k2;
+ k2cpt++;
+ /* printf("iind=%ld dateintmean = %lf dateintsum=%lf k2cpt=%lf k2=%lf\n",iind, dateintsum/k2cpt, dateintsum,k2cpt, k2); */
+ }
+ } /* end bool 2 */
+ } /* end m */
} /* end bool */
} /* end iind = 1 to imx */
/* prop[s][age] is feeded for any initial and valid live state as well as
freq[s1][s2][age] at single age of beginning the transition, for a combination j1 */
-
-
+
+
/* fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
pstamp(ficresp);
/* if (ncoveff>0) { */
fprintf(ficresphtm, "\n<br/><br/><h3>********** Variable ");
fprintf(ficresphtmfr, "\n<br/><br/><h3>********** Variable ");
for (z1=1; z1<=cptcoveff; z1++){
- fprintf(ficresp, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
- fprintf(ficresphtm, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
- fprintf(ficresphtmfr, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresp, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresphtm, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresphtmfr, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
}
fprintf(ficresp, "**********\n#");
fprintf(ficresphtm, "**********</h3>\n");
}
fprintf(ficresp, "\n");
fprintf(ficresphtm, "\n");
-
+
/* Header of frequency table by age */
fprintf(ficresphtmfr,"<table style=\"text-align:center; border: 1px solid\">");
fprintf(ficresphtmfr,"<th>Age</th> ");
for(jk=-1; jk <=nlstate+ndeath; jk++){
for(m=-1; m <=nlstate+ndeath; m++){
- if(jk!=0 && m!=0)
- fprintf(ficresphtmfr,"<th>%d%d</th> ",jk,m);
+ if(jk!=0 && m!=0)
+ fprintf(ficresphtmfr,"<th>%d%d</th> ",jk,m);
}
}
fprintf(ficresphtmfr, "\n");
-
+
/* For each age */
for(iage=iagemin; iage <= iagemax+3; iage++){
fprintf(ficresphtm,"<tr>");
if(iage==iagemax+1){
- fprintf(ficlog,"1");
- fprintf(ficresphtmfr,"<tr><th>0</th> ");
+ fprintf(ficlog,"1");
+ fprintf(ficresphtmfr,"<tr><th>0</th> ");
}else if(iage==iagemax+2){
- fprintf(ficlog,"0");
- fprintf(ficresphtmfr,"<tr><th>Unknown</th> ");
+ fprintf(ficlog,"0");
+ fprintf(ficresphtmfr,"<tr><th>Unknown</th> ");
}else if(iage==iagemax+3){
- fprintf(ficlog,"Total");
- fprintf(ficresphtmfr,"<tr><th>Total</th> ");
+ fprintf(ficlog,"Total");
+ fprintf(ficresphtmfr,"<tr><th>Total</th> ");
}else{
- if(first==1){
- first=0;
- printf("See log file for details...\n");
- }
- fprintf(ficresphtmfr,"<tr><th>%d</th> ",iage);
- fprintf(ficlog,"Age %d", iage);
+ if(first==1){
+ first=0;
+ printf("See log file for details...\n");
+ }
+ fprintf(ficresphtmfr,"<tr><th>%d</th> ",iage);
+ fprintf(ficlog,"Age %d", iage);
}
for(jk=1; jk <=nlstate ; jk++){
- for(m=-1, pp[jk]=0; m <=nlstate+ndeath ; m++)
- pp[jk] += freq[jk][m][iage];
+ for(m=-1, pp[jk]=0; m <=nlstate+ndeath ; m++)
+ pp[jk] += freq[jk][m][iage];
}
for(jk=1; jk <=nlstate ; jk++){
- for(m=-1, pos=0; m <=0 ; m++)
- pos += freq[jk][m][iage];
- if(pp[jk]>=1.e-10){
- if(first==1){
- printf(" %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
- }
- fprintf(ficlog," %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
- }else{
- if(first==1)
- printf(" %d.=%.0f loss[%d]=NaNQ%%",jk,pp[jk],jk);
- fprintf(ficlog," %d.=%.0f loss[%d]=NaNQ%%",jk,pp[jk],jk);
- }
+ for(m=-1, pos=0; m <=0 ; m++)
+ pos += freq[jk][m][iage];
+ if(pp[jk]>=1.e-10){
+ if(first==1){
+ printf(" %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
+ }
+ fprintf(ficlog," %d.=%.0f loss[%d]=%.1f%%",jk,pp[jk],jk,100*pos/pp[jk]);
+ }else{
+ if(first==1)
+ printf(" %d.=%.0f loss[%d]=NaNQ%%",jk,pp[jk],jk);
+ fprintf(ficlog," %d.=%.0f loss[%d]=NaNQ%%",jk,pp[jk],jk);
+ }
}
-
+
for(jk=1; jk <=nlstate ; jk++){
- /* posprop[jk]=0; */
- for(m=0, pp[jk]=0; m <=nlstate+ndeath; m++)/* Summing on all ages */
- pp[jk] += freq[jk][m][iage];
+ /* posprop[jk]=0; */
+ for(m=0, pp[jk]=0; m <=nlstate+ndeath; m++)/* Summing on all ages */
+ pp[jk] += freq[jk][m][iage];
} /* pp[jk] is the total number of transitions starting from state jk and any ending status until this age */
-
+
for(jk=1,pos=0, pospropta=0.; jk <=nlstate ; jk++){
- pos += pp[jk]; /* pos is the total number of transitions until this age */
- posprop[jk] += prop[jk][iage]; /* prop is the number of transitions from a live state
- from jk at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */
- pospropta += prop[jk][iage]; /* prop is the number of transitions from a live state
- from jk at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */
+ pos += pp[jk]; /* pos is the total number of transitions until this age */
+ posprop[jk] += prop[jk][iage]; /* prop is the number of transitions from a live state
+ from jk at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */
+ pospropta += prop[jk][iage]; /* prop is the number of transitions from a live state
+ from jk at age iage prop[s[m][iind]][(int)agev[m][iind]] += weight[iind] */
}
for(jk=1; jk <=nlstate ; jk++){
- if(pos>=1.e-5){
- if(first==1)
- printf(" %d.=%.0f prev[%d]=%.1f%%",jk,pp[jk],jk,100*pp[jk]/pos);
- fprintf(ficlog," %d.=%.0f prev[%d]=%.1f%%",jk,pp[jk],jk,100*pp[jk]/pos);
- }else{
- if(first==1)
- printf(" %d.=%.0f prev[%d]=NaNQ%%",jk,pp[jk],jk);
- fprintf(ficlog," %d.=%.0f prev[%d]=NaNQ%%",jk,pp[jk],jk);
- }
- if( iage <= iagemax){
- if(pos>=1.e-5){
- fprintf(ficresp," %d %.5f %.0f %.0f",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta);
- fprintf(ficresphtm,"<th>%d</th><td>%.5f</td><td>%.0f</td><td>%.0f</td>",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta);
- /*probs[iage][jk][j1]= pp[jk]/pos;*/
- /*printf("\niage=%d jk=%d j1=%d %.5f %.0f %.0f %f",iage,jk,j1,pp[jk]/pos, pp[jk],pos,probs[iage][jk][j1]);*/
- }
- else{
- fprintf(ficresp," %d NaNq %.0f %.0f",iage,prop[jk][iage],pospropta);
- fprintf(ficresphtm,"<th>%d</th><td>NaNq</td><td>%.0f</td><td>%.0f</td>",iage, prop[jk][iage],pospropta);
- }
- }
- pospropt[jk] +=posprop[jk];
+ if(pos>=1.e-5){
+ if(first==1)
+ printf(" %d.=%.0f prev[%d]=%.1f%%",jk,pp[jk],jk,100*pp[jk]/pos);
+ fprintf(ficlog," %d.=%.0f prev[%d]=%.1f%%",jk,pp[jk],jk,100*pp[jk]/pos);
+ }else{
+ if(first==1)
+ printf(" %d.=%.0f prev[%d]=NaNQ%%",jk,pp[jk],jk);
+ fprintf(ficlog," %d.=%.0f prev[%d]=NaNQ%%",jk,pp[jk],jk);
+ }
+ if( iage <= iagemax){
+ if(pos>=1.e-5){
+ fprintf(ficresp," %d %.5f %.0f %.0f",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta);
+ fprintf(ficresphtm,"<th>%d</th><td>%.5f</td><td>%.0f</td><td>%.0f</td>",iage,prop[jk][iage]/pospropta, prop[jk][iage],pospropta);
+ /*probs[iage][jk][j1]= pp[jk]/pos;*/
+ /*printf("\niage=%d jk=%d j1=%d %.5f %.0f %.0f %f",iage,jk,j1,pp[jk]/pos, pp[jk],pos,probs[iage][jk][j1]);*/
+ }
+ else{
+ fprintf(ficresp," %d NaNq %.0f %.0f",iage,prop[jk][iage],pospropta);
+ fprintf(ficresphtm,"<th>%d</th><td>NaNq</td><td>%.0f</td><td>%.0f</td>",iage, prop[jk][iage],pospropta);
+ }
+ }
+ pospropt[jk] +=posprop[jk];
} /* end loop jk */
/* pospropt=0.; */
for(jk=-1; jk <=nlstate+ndeath; jk++){
- for(m=-1; m <=nlstate+ndeath; m++){
- if(freq[jk][m][iage] !=0 ) { /* minimizing output */
- if(first==1){
- printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]);
- }
- fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iage]);
- }
- if(jk!=0 && m!=0)
- fprintf(ficresphtmfr,"<td>%.0f</td> ",freq[jk][m][iage]);
- }
+ for(m=-1; m <=nlstate+ndeath; m++){
+ if(freq[jk][m][iage] !=0 ) { /* minimizing output */
+ if(first==1){
+ printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]);
+ }
+ fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iage]);
+ }
+ if(jk!=0 && m!=0)
+ fprintf(ficresphtmfr,"<td>%.0f</td> ",freq[jk][m][iage]);
+ }
} /* end loop jk */
posproptt=0.;
for(jk=1; jk <=nlstate; jk++){
- posproptt += pospropt[jk];
+ posproptt += pospropt[jk];
}
fprintf(ficresphtmfr,"</tr>\n ");
if(iage <= iagemax){
- fprintf(ficresp,"\n");
- fprintf(ficresphtm,"</tr>\n");
+ fprintf(ficresp,"\n");
+ fprintf(ficresphtm,"</tr>\n");
}
if(first==1)
- printf("Others in log...\n");
+ printf("Others in log...\n");
fprintf(ficlog,"\n");
} /* end loop age iage */
fprintf(ficresphtm,"<tr><th>Tot</th>");
for(jk=1; jk <=nlstate ; jk++){
if(posproptt < 1.e-5){
- fprintf(ficresphtm,"<td>Nanq</td><td>%.0f</td><td>%.0f</td>",pospropt[jk],posproptt);
+ fprintf(ficresphtm,"<td>Nanq</td><td>%.0f</td><td>%.0f</td>",pospropt[jk],posproptt);
}else{
- fprintf(ficresphtm,"<td>%.5f</td><td>%.0f</td><td>%.0f</td>",pospropt[jk]/posproptt,pospropt[jk],posproptt);
+ fprintf(ficresphtm,"<td>%.5f</td><td>%.0f</td><td>%.0f</td>",pospropt[jk]/posproptt,pospropt[jk],posproptt);
}
}
fprintf(ficresphtm,"</tr>\n");
fprintf(ficresphtmfr,"</table>\n");
} /* end selected combination of covariate j1 */
dateintmean=dateintsum/k2cpt;
-
+
fclose(ficresp);
fclose(ficresphtm);
fclose(ficresphtmfr);
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 */
-
- /* 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;
+ 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 */
+
+ /* 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;
+ ij=0;
+ for(i=0; i<=1;i++){
+ nbcode[Tvar[k]][++ij]=i;
+ }
+ break;
default:
- break;
+ break;
} /* end switch */
} /* end dummy test */
++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[k]=Tvar[k]- ncovcol-nqv; /* Inverse TmodelInvind[2=V4]=2 second dummy varying cov (V4)4-1-1 {0, 2, 1, } TmodelInvind[3]=1 */
+ 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 *\/ */
+ /* }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?*/
+ * 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;
if ((strpt=strstr(model,"age*age")) !=0){
printf(" strpt=%s, model=%s\n",strpt, model);
if(strpt != model){
- printf("Error in model: 'model=%s'; 'age*age' should in first place before other covariates\n \
+ printf("Error in model: 'model=%s'; 'age*age' should in first place before other covariates\n \
'model=1+age+age*age+V1.' or 'model=1+age+age*age+V1+V1*age.', please swap as well as \n \
corresponding column of parameters.\n",model);
- fprintf(ficlog,"Error in model: 'model=%s'; 'age*age' should in first place before other covariates\n \
+ fprintf(ficlog,"Error in model: 'model=%s'; 'age*age' should in first place before other covariates\n \
'model=1+age+age*age+V1.' or 'model=1+age+age*age+V1+V1*age.', please swap as well as \n \
corresponding column of parameters.\n",model); fflush(ficlog);
- return 1;
+ return 1;
}
nagesqr=1;
if (strstr(model,"+age*age") !=0)
- substrchaine(modelsav, model, "+age*age");
+ substrchaine(modelsav, model, "+age*age");
else if (strstr(model,"age*age+") !=0)
- substrchaine(modelsav, model, "age*age+");
+ substrchaine(modelsav, model, "age*age+");
else
- substrchaine(modelsav, model, "age*age");
+ substrchaine(modelsav, model, "age*age");
}else
nagesqr=0;
if (strlen(modelsav) >1){
Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product\n",model);
for(k=1, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0;k<=cptcovt; k++){ /* or cptocvt */
- if (Tvar[k] <=ncovcol && (Typevar[k]==0 || Typevar[k]==2)){ /* Simple or product fixed dummy covariatee */
+ if (Tvar[k] <=ncovcol && (Typevar[k]==0 || Typevar[k]==2)){ /* Simple or product fixed dummy (<=ncovcol) covariates */
Fixed[k]= 0;
Dummy[k]= 0;
ncoveff++;
+ modell[k].maintype= FTYPE;
TvarFD[ncoveff]=Tvar[k]; /* TvarFD[1]=V1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
- TvarFDind[ncoveff]=Tvar[k]; /* TvarFDind[1]=9 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
- }else if( Tvar[k] <=ncovcol+nqv && Typevar[k]==0){ /* Remind that product Vn*Vm are added in k*/
+ TvarFDind[ncoveff]=k; /* TvarFDind[1]=9 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */
+ }else if( Tvar[k] <=ncovcol+nqv && Typevar[k]==0){ /* Remind that product Vn*Vm are added in k*/ /* Only simple fixed quantitative variable */
Fixed[k]= 0;
Dummy[k]= 1;
nqfveff++;
- TvarFQ[nqfveff]=Tvar[k]; /* TvarFQ[1]=V2 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
+ modell[k].maintype= FTYPE;
+ modell[k].subtype= FQ;
+ 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){
Fixed[k]= 1;
Dummy[k]= 0;
ntveff++; /* Only simple time varying dummy variable */
- TvarVD[ntvveff]=Tvar[k]; /* TvarVD[1]=V5 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
- TvarVDind[ntveff++]=k; /* TvarVDind[1]=1 in V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1 */ /* Only simple fixed quantitative variable */
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VD;
+ 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);
printf("Quasi TmodelInvind[%d]=%d\n",k,Tvar[k]- ncovcol-nqv);
- }else if( Tvar[k] <=ncovcol+nqv+ntv+nqtv && Typevar[k]==0){
- Fixed[k]= 1;
- Dummy[k]= 1;
- TmodelInvQind[++nqtveff]=Tvar[k]- ncovcol-nqv-ntv;/* Only simple time varying quantitative variable */
- /* Tmodeliqind[k]=nqtveff;/\* Only simple time varying quantitative variable *\/ */
- printf("Quasi TmodelQind[%d]=%d,Tvar[TmodelQind[%d]]=V%d, ncovcol=%d, nqv=%d, ntv=%d,Tvar[k]- ncovcol-nqv-ntv=%d\n",nqtveff,k,nqtveff,Tvar[k], ncovcol, nqv, ntv, Tvar[k]- ncovcol-nqv-ntv);
+ }else if( Tvar[k] <=ncovcol+nqv+ntv+nqtv && Typevar[k]==0){ /* Only simple time varying quantitative variable V5*/
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ nqtveff++;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VQ;
+ 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 */
+ /* Tmodeliqind[k]=nqtveff;/\* Only simple time varying quantitative variable *\/ */
+ printf("Quasi TmodelQind[%d]=%d,Tvar[TmodelQind[%d]]=V%d, ncovcol=%d, nqv=%d, ntv=%d,Tvar[k]- ncovcol-nqv-ntv=%d\n",nqtveff,k,nqtveff,Tvar[k], ncovcol, nqv, ntv, Tvar[k]- ncovcol-nqv-ntv);
printf("Quasi TmodelInvQind[%d]=%d\n",k,Tvar[k]- ncovcol-nqv-ntv);
}else if (Typevar[k] == 1) { /* product with age */
- if (Tvar[k] <=ncovcol ){ /* Simple or product fixed dummy covariatee */
- Fixed[k]= 2;
- Dummy[k]= 2;
- /* ncoveff++; */
+ if (Tvar[k] <=ncovcol ){ /* Product age with fixed dummy covariatee */
+ Fixed[k]= 2;
+ Dummy[k]= 2;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APFD;
+ /* ncoveff++; */
}else if( Tvar[k] <=ncovcol+nqv) { /* Remind that product Vn*Vm are added in k*/
- Fixed[k]= 2;
- Dummy[k]= 3;
- /* nqfveff++; /\* Only simple fixed quantitative variable *\/ */
+ Fixed[k]= 2;
+ Dummy[k]= 3;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APFQ; /* Product age * fixed quantitative */
+ /* nqfveff++; /\* Only simple fixed quantitative variable *\/ */
}else if( Tvar[k] <=ncovcol+nqv+ntv ){
- Fixed[k]= 3;
- Dummy[k]= 2;
- /* ntveff++; /\* Only simple time varying dummy variable *\/ */
+ Fixed[k]= 3;
+ Dummy[k]= 2;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APVD; /* Product age * varying dummy */
+ /* ntveff++; /\* Only simple time varying dummy variable *\/ */
}else if( Tvar[k] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 3;
- Dummy[k]= 3;
- /* nqtveff++;/\* Only simple time varying quantitative variable *\/ */
+ Fixed[k]= 3;
+ Dummy[k]= 3;
+ modell[k].maintype= ATYPE;
+ modell[k].subtype= APVQ; /* Product age * varying quantitative */
+ /* nqtveff++;/\* Only simple time varying quantitative variable *\/ */
}
}else if (Typevar[k] == 2) { /* product without age */
k1=Tposprod[k];
if(Tvard[k1][1] <=ncovcol){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 1;
- Dummy[k]= 0;
- }else if(Tvard[k1][2] <=ncovcol+nqv){
- Fixed[k]= 0; /* or 2 ?*/
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 0;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 1;
+ Dummy[k]= 0;
+ modell[k].maintype= FTYPE;
+ modell[k].subtype= FPDD; /* Product fixed dummy * fixed dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv){
+ Fixed[k]= 0; /* or 2 ?*/
+ Dummy[k]= 1;
+ modell[k].maintype= FTYPE;
+ modell[k].subtype= FPDQ; /* Product fixed dummy * fixed quantitative */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 0;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDD; /* Product fixed dummy * varying dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product fixed dummy * varying quantitative */
+ }
}else if(Tvard[k1][1] <=ncovcol+nqv){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 0; /* or 2 ?*/
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv){
- Fixed[k]= 0; /* or 2 ?*/
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 0; /* or 2 ?*/
+ Dummy[k]= 1;
+ modell[k].maintype= FTYPE;
+ modell[k].subtype= FPDQ; /* Product fixed quantitative * fixed dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product fixed quantitative * varying dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPQQ; /* Product fixed quantitative * varying quantitative */
+ }
}else if(Tvard[k1][1] <=ncovcol+nqv+ntv){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 0;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDD; /* Product time varying dummy * fixed dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying dummy * fixed quantitative */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 0;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDD; /* Product time varying dummy * time varying dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying dummy * time varying quantitative */
+ }
}else if(Tvard[k1][1] <=ncovcol+nqv+ntv+nqtv){
- if(Tvard[k1][2] <=ncovcol){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
- Fixed[k]= 1;
- Dummy[k]= 1;
- }
+ if(Tvard[k1][2] <=ncovcol){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying quantitative * fixed dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPQQ; /* Product time varying quantitative * fixed quantitative */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPDQ; /* Product time varying quantitative * time varying dummy */
+ }else if(Tvard[k1][2] <=ncovcol+nqv+ntv+nqtv){
+ Fixed[k]= 1;
+ Dummy[k]= 1;
+ modell[k].maintype= VTYPE;
+ modell[k].subtype= VPQQ; /* Product time varying quantitative * time varying quantitative */
+ }
}else{
- printf("Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
- fprintf(ficlog,"Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
+ printf("Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
+ fprintf(ficlog,"Error unknown type of covariate: Tvard[%d][1]=%d,Tvard[%d][2]=%d\n",k1,Tvard[k1][1],k1,Tvard[k1][2]);
} /* end k1 */
}else{
printf("Error, current version can't treat for performance reasons, Tvar[%d]=%d, Typevar[%d]=%d\n", k, Tvar[k], k, Typevar[k]);
fprintf(ficlog,"Error, current version can't treat for performance reasons, Tvar[%d]=%d, Typevar[%d]=%d\n", k, Tvar[k], k, Typevar[k]);
}
printf("Decodemodel, k=%d, Tvar[%d]=V%d,Typevar=%d, Fixed=%d, Dummy=%d\n",k, k,Tvar[k],Typevar[k],Fixed[k],Dummy[k]);
+ printf(" modell[%d].maintype=%d, modell[%d].subtype=%d\n",k,modell[k].maintype,k,modell[k].subtype);
fprintf(ficlog,"Decodemodel, k=%d, Tvar[%d]=V%d,Typevar=%d, Fixed=%d, Dummy=%d\n",k, k,Tvar[k],Typevar[k],Fixed[k],Dummy[k]);
}
/* Searching for doublons in the model */
for(k1=1; k1<= cptcovt;k1++){
for(k2=1; k2 <k1;k2++){
if((Typevar[k1]==Typevar[k2]) && (Fixed[Tvar[k1]]==Fixed[Tvar[k2]]) && (Dummy[Tvar[k1]]==Dummy[Tvar[k2]] )){
- if((Typevar[k1] == 0 || Typevar[k1] == 1)){ /* Simple or age product */
- if(Tvar[k1]==Tvar[k2]){
- printf("Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]);
- fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); fflush(ficlog);
- return(1);
- }
- }else if (Typevar[k1] ==2){
- k3=Tposprod[k1];
- k4=Tposprod[k2];
- if( ((Tvard[k3][1]== Tvard[k4][1])&&(Tvard[k3][2]== Tvard[k4][2])) || ((Tvard[k3][1]== Tvard[k4][2])&&(Tvard[k3][2]== Tvard[k4][1])) ){
- printf("Error duplication in the model=%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]);
- fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); fflush(ficlog);
- return(1);
- }
- }
+ if((Typevar[k1] == 0 || Typevar[k1] == 1)){ /* Simple or age product */
+ if(Tvar[k1]==Tvar[k2]){
+ printf("Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]);
+ fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); fflush(ficlog);
+ return(1);
+ }
+ }else if (Typevar[k1] ==2){
+ k3=Tposprod[k1];
+ k4=Tposprod[k2];
+ if( ((Tvard[k3][1]== Tvard[k4][1])&&(Tvard[k3][2]== Tvard[k4][2])) || ((Tvard[k3][1]== Tvard[k4][2])&&(Tvard[k3][2]== Tvard[k4][1])) ){
+ printf("Error duplication in the model=%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]);
+ fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, V%d*V%d, Typevar=%d, Fixed=%d, Dummy=%d\n",model, k1,k2, Tvard[k3][1], Tvard[k3][2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); fflush(ficlog);
+ return(1);
+ }
+ }
}
}
}
k=2 V1 Tvar[k=2]= 1 (from V1)
k=1 Tvar[1]=2 (from V2)
*/
- Tvar=ivector(1,NCOVMAX); /* Was 15 changed to NCOVMAX. */
- Tvarsel=ivector(1,NCOVMAX); /* */
+
+ Tvar=ivector(1,NCOVMAX); /* Was 15 changed to NCOVMAX. */
+ TvarFD=ivector(1,NCOVMAX); /* */
+ TvarFDind=ivector(1,NCOVMAX); /* */
+ TvarFQ=ivector(1,NCOVMAX); /* */
+ TvarFQind=ivector(1,NCOVMAX); /* */
+ TvarVD=ivector(1,NCOVMAX); /* */
+ TvarVDind=ivector(1,NCOVMAX); /* */
+ TvarVQ=ivector(1,NCOVMAX); /* */
+ TvarVQind=ivector(1,NCOVMAX); /* */
+
Tvalsel=vector(1,NCOVMAX); /* */
Typevar=ivector(-1,NCOVMAX); /* -1 to 2 */
Fixed=ivector(-1,NCOVMAX); /* -1 to 3 */
free_ivector(Fixed,-1,NCOVMAX);
free_ivector(Typevar,-1,NCOVMAX);
free_ivector(Tvar,1,NCOVMAX);
+ free_ivector(TvarFD,1,NCOVMAX);
+ free_ivector(TvarFDind,1,NCOVMAX);
+ free_ivector(TvarFQ,1,NCOVMAX);
+ free_ivector(TvarFQind,1,NCOVMAX);
+ free_ivector(TvarVD,1,NCOVMAX);
+ free_ivector(TvarVDind,1,NCOVMAX);
+ free_ivector(TvarVQ,1,NCOVMAX);
+ free_ivector(TvarVQind,1,NCOVMAX);
free_ivector(Tvarsel,1,NCOVMAX);
free_vector(Tvalsel,1,NCOVMAX);
free_ivector(Tposprod,1,NCOVMAX);
git://henry.ined.fr / .git / commitdiff