\n- Convergence to period (stable) prevalence in state %d. Or probability for a person being in state (1 to %d) at different ages, to be in state %d some years after. %s_%d-%d-%d.svg
\
![](\"%s_%d-%d-%d.svg\")
", cpt, nlstate, cpt, subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres);
}
- if(backcast==1){
+ if(prevbcast==1){
/* Backward prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
fprintf(fichtm,"\n- Convergence to mixed (stable) back prevalence in state %d. Or probability for a person to be in state %d at a younger age, knowing that she/he was in state (1 to %d) at different older ages. %s_%d-%d-%d.svg
\ @@ -6964,24 +7094,28 @@ divided by h: hPij if(prevfcast==1){ /* Projection of prevalence up to period (forward stable) prevalence in each health state */ for(cpt=1; cpt<=nlstate;cpt++){ - fprintf(fichtm,"
\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f and mobil_average=%d), from year %.1f up to year %.1f tending to period (stable) forward prevalence in state %d. Or probability to be in state %d being in an observed weighted state (from 1 to %d). %s_%d-%d-%d.svg
\ -
", dateprev1, dateprev2, mobilavproj, dateproj1, dateproj2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres);
+ fprintf(fichtm,"\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f and mobil_average=%d), from year %.1f up to year %.1f tending to period (stable) forward prevalence in state %d. Or probability to be in state %d being in an observed weighted state (from 1 to %d). %s_%d-%d-%d.svg", dateprev1, dateprev2, mobilavproj, dateprojd, dateprojf, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres); + fprintf(fichtm," (data from text file %s.txt)\n
",subdirf2(optionfilefiname,"F_"),subdirf2(optionfilefiname,"F_")); + fprintf(fichtm,"
",
+ subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres);
}
}
- if(backcast==1){
+ if(prevbcast==1){
/* Back projection of prevalence up to stable (mixed) back-prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
fprintf(fichtm,"\n- Back projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f and mobil_average=%d), \ from year %.1f up to year %.1f (probably close to stable [mixed] back prevalence in state %d (randomness in cross-sectional prevalence is not taken into \ account but can visually be appreciated). Or probability to have been in an state %d, knowing that the person was in either state (1 or %d) \ -with weights corresponding to observed prevalence at different ages. %s_%d-%d-%d.svg
\ -
", dateprev1, dateprev2, mobilavproj, dateback1, dateback2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJB_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJB_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJB_"),cpt,k1,nres);
+with weights corresponding to observed prevalence at different ages. %s_%d-%d-%d.svg", dateprev1, dateprev2, mobilavproj, dateback1, dateback2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJB_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJB_"),cpt,k1,nres);
+ fprintf(fichtm," (data from text file %s.txt)\n",subdirf2(optionfilefiname,"FB_"),subdirf2(optionfilefiname,"FB_")); + fprintf(fichtm,"
", subdirf2(optionfilefiname,"PROJB_"),cpt,k1,nres);
}
}
for(cpt=1; cpt<=nlstate;cpt++) {
- fprintf(fichtm,"\n- Life expectancy by health state (%d) at initial age and its decomposition into health expectancies in each alive state (1 to %d) (or area under each survival functions): %s_%d-%d-%d.svg
\ -
",cpt,nlstate,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres);
+ fprintf(fichtm,"\n- Life expectancy by health state (%d) at initial age and its decomposition into health expectancies in each alive state (1 to %d) (or area under each survival functions): %s_%d-%d-%d.svg",cpt,nlstate,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres,subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres); + fprintf(fichtm," (data from text file %s.txt)\n
",subdirf2(optionfilefiname,"E_"),subdirf2(optionfilefiname,"E_")); + fprintf(fichtm,"
", subdirf2(optionfilefiname,"EXP_"),cpt,k1,nres );
}
/* } /\* end i1 *\/ */
}/* End k1 */
@@ -7064,15 +7198,17 @@ See page 'Matrix of variance-covariance
}
for(cpt=1; cpt<=nlstate;cpt++) {
fprintf(fichtm,"\n- Observed (cross-sectional with mov_average=%d) and period (incidence based) \ -prevalence (with 95%% confidence interval) in state (%d): %s_%d-%d-%d.svg\n
\ -
",mobilav,cpt,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres);
+prevalence (with 95%% confidence interval) in state (%d): %s_%d-%d-%d.svg",mobilav,cpt,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres);
+ fprintf(fichtm," (data from text file %s)\n ",subdirf2(fileresu,"VPL_"),subdirf2(fileresu,"VPL_")); + fprintf(fichtm,"
",subdirf2(optionfilefiname,"V_"), cpt,k1,nres);
}
fprintf(fichtm,"\n- Total life expectancy by age and \ -health expectancies in states (1) and (2). If popbased=1 the smooth (due to the model) \ +health expectancies in each live states (1 to %d). If popbased=1 the smooth (due to the model) \ true period expectancies (those weighted with period prevalences are also\ drawn in addition to the population based expectancies computed using\ - observed and cahotic prevalences: %s_%d-%d.svg\n
\ -
![](\"%s_%d-%d.svg\")
",subdirf2(optionfilefiname,"E_"),k1,nres,subdirf2(optionfilefiname,"E_"),k1,nres,subdirf2(optionfilefiname,"E_"),k1,nres);
+ observed and cahotic prevalences: %s_%d-%d.svg",nlstate, subdirf2(optionfilefiname,"E_"),k1,nres,subdirf2(optionfilefiname,"E_"),k1,nres);
+ fprintf(fichtm," (data from text file %s.txt) \n",subdirf2(optionfilefiname,"T_"),subdirf2(optionfilefiname,"T_")); + fprintf(fichtm,"
",subdirf2(optionfilefiname,"E_"),k1,nres);
/* } /\* end i1 *\/ */
}/* End k1 */
}/* End nres */
@@ -7081,7 +7217,7 @@ true period expectancies (those weighted
}
/******************* Gnuplot file **************/
-void printinggnuplot(char fileresu[], char optionfilefiname[], double ageminpar, double agemaxpar, double bage, double fage , int prevfcast, int backcast, char pathc[], double p[], int offyear, int offbyear){
+void printinggnuplot(char fileresu[], char optionfilefiname[], double ageminpar, double agemaxpar, double bage, double fage , int prevfcast, int prevbcast, char pathc[], double p[], int offyear, int offbyear){
char dirfileres[132],optfileres[132];
char gplotcondition[132], gplotlabel[132];
@@ -7235,7 +7371,7 @@ void printinggnuplot(char fileresu[], ch
} /* end covariate */
} /* end if no covariate */
- if(backcast==1){ /* We need to get the corresponding values of the covariates involved in this combination k1 */
+ if(prevbcast==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 1:((",subdirf2(fileresu,"PLB_")); /* Age is in 1, nres in 2 to be fixed */
if(cptcoveff ==0){
@@ -7262,7 +7398,7 @@ void printinggnuplot(char fileresu[], ch
}
} /* end covariate */
} /* end if no covariate */
- if(backcast == 1){
+ if(prevbcast == 1){
fprintf(ficgp,", \"%s\" every :::%d::%d u 1:($2==%d ? $3:1/0) \"%%lf %%lf",subdirf2(fileresu,"VBL_"),nres-1,nres-1,nres);
/* k1-1 error should be nres-1*/
for (i=1; i<= nlstate ; i ++) {
@@ -7281,7 +7417,7 @@ void printinggnuplot(char fileresu[], ch
}
fprintf(ficgp,"\" t\"\" w l lt 4");
} /* end if backprojcast */
- } /* end if backcast */
+ } /* end if prevbcast */
/* fprintf(ficgp,"\nset out ;unset label;\n"); */
fprintf(ficgp,"\nset out ;unset title;\n");
} /* nres */
@@ -7569,7 +7705,7 @@ set ter svg size 640, 480\nunset log y\n
/* 7eme */
- if(backcast == 1){
+ if(prevbcast == 1){
/* CV backward prevalence for each covariate */
for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
@@ -7621,7 +7757,7 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficgp,"\nset out; unset label;\n");
} /* end cpt state*/
} /* end covariate */
- } /* End if backcast */
+ } /* End if prevbcast */
/* 8eme */
if(prevfcast==1){
@@ -7737,7 +7873,7 @@ set ter svg size 640, 480\nunset log y\n
} /* end covariate */
} /* End if prevfcast */
- if(backcast==1){
+ if(prevbcast==1){
/* Back projection from cross-sectional to stable (mixed) for each covariate */
for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
@@ -7850,7 +7986,7 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficgp,"\nset out; unset label;\n");
} /* end cpt state*/
} /* end covariate */
- } /* End if backcast */
+ } /* End if prevbcast */
/* 9eme writing MLE parameters */
@@ -8298,16 +8434,21 @@ set ter svg size 640, 480\nunset log y\n
}/* End movingaverage */
+
/************** Forecasting ******************/
- void prevforecast(char fileres[], double anproj1, double mproj1, double jproj1, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double ***prev, double bage, double fage, int firstpass, int lastpass, double anproj2, double p[], int cptcoveff){
- /* proj1, year, month, day of starting projection
+/* void prevforecast(char fileres[], double dateintmean, double anprojd, double mprojd, double jprojd, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double ***prev, double bage, double fage, int firstpass, int lastpass, double anprojf, double p[], int cptcoveff)*/
+void prevforecast(char fileres[], double dateintmean, double dateprojd, double dateprojf, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double ***prev, double bage, double fage, int firstpass, int lastpass, double p[], int cptcoveff){
+ /* dateintemean, mean date of interviews
+ dateprojd, year, month, day of starting projection
+ dateprojf date of end of projection;year of end of projection (same day and month as proj1).
agemin, agemax range of age
dateprev1 dateprev2 range of dates during which prevalence is computed
- anproj2 year of en of projection (same day and month as proj1).
*/
+ /* double anprojd, mprojd, jprojd; */
+ /* double anprojf, mprojf, jprojf; */
int yearp, stepsize, hstepm, nhstepm, j, k, cptcod, i, h, i1, k4, nres=0;
double agec; /* generic age */
- double agelim, ppij, yp,yp1,yp2,jprojmean,mprojmean,anprojmean;
+ double agelim, ppij, yp,yp1,yp2;
double *popeffectif,*popcount;
double ***p3mat;
/* double ***mobaverage; */
@@ -8344,22 +8485,27 @@ set ter svg size 640, 480\nunset log y\n
if(estepm > stepm){ /* Yes every two year */
stepsize=2;
}
+ hstepm=hstepm/stepm;
+
+
+ /* yp1=modf(dateintmean,&yp);/\* extracts integral of datemean in yp and */
+ /* fractional in yp1 *\/ */
+ /* aintmean=yp; */
+ /* yp2=modf((yp1*12),&yp); */
+ /* mintmean=yp; */
+ /* yp1=modf((yp2*30.5),&yp); */
+ /* jintmean=yp; */
+ /* if(jintmean==0) jintmean=1; */
+ /* if(mintmean==0) mintmean=1; */
- hstepm=hstepm/stepm;
- yp1=modf(dateintmean,&yp);/* extracts integral of datemean in yp and
- fractional in yp1 */
- anprojmean=yp;
- yp2=modf((yp1*12),&yp);
- mprojmean=yp;
- yp1=modf((yp2*30.5),&yp);
- jprojmean=yp;
- if(jprojmean==0) jprojmean=1;
- if(mprojmean==0) jprojmean=1;
+ /* date2dmy(dateintmean,&jintmean,&mintmean,&aintmean); */
+ /* date2dmy(dateprojd,&jprojd, &mprojd, &anprojd); */
+ /* date2dmy(dateprojf,&jprojf, &mprojf, &anprojf); */
i1=pow(2,cptcoveff);
if (cptcovn < 1){i1=1;}
- fprintf(ficresf,"# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jprojmean,mprojmean,anprojmean,dateintmean,dateprev1,dateprev2);
+ fprintf(ficresf,"# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jintmean,mintmean,aintmean,dateintmean,dateprev1,dateprev2);
fprintf(ficresf,"#****** Routine prevforecast **\n");
@@ -8385,9 +8531,9 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficresf," p%d%d",i,j);
fprintf(ficresf," wp.%d",j);
}
- for (yearp=0; yearp<=(anproj2-anproj1);yearp +=stepsize) {
+ for (yearp=0; yearp<=(anprojf-anprojd);yearp +=stepsize) {
fprintf(ficresf,"\n");
- fprintf(ficresf,"\n# Forecasting at date %.lf/%.lf/%.lf ",jproj1,mproj1,anproj1+yearp);
+ fprintf(ficresf,"\n# Forecasting at date %.lf/%.lf/%.lf ",jprojd,mprojd,anprojd+yearp);
/* for (agec=fage; agec>=(ageminpar-1); agec--){ */
for (agec=fage; agec>=(bage); agec--){
nhstepm=(int) rint((agelim-agec)*YEARM/stepm);
@@ -8405,7 +8551,7 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficresf,"\n");
for(j=1;j<=cptcoveff;j++)
fprintf(ficresf,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- fprintf(ficresf,"%.f %.f ",anproj1+yearp,agec+h*hstepm/YEARM*stepm);
+ fprintf(ficresf,"%.f %.f ",anprojd+yearp,agec+h*hstepm/YEARM*stepm);
for(j=1; j<=nlstate+ndeath;j++) {
ppij=0.;
@@ -8433,8 +8579,9 @@ set ter svg size 640, 480\nunset log y\n
}
/************** Back Forecasting ******************/
- void prevbackforecast(char fileres[], double ***prevacurrent, double anback1, double mback1, double jback1, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double bage, double fage, int firstpass, int lastpass, double anback2, double p[], int cptcoveff){
- /* back1, year, month, day of starting backection
+ /* void prevbackforecast(char fileres[], double ***prevacurrent, double anback1, double mback1, double jback1, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double bage, double fage, int firstpass, int lastpass, double anback2, double p[], int cptcoveff){ */
+ void prevbackforecast(char fileres[], double ***prevacurrent, double dateintmean, double dateprojd, double dateprojf, double ageminpar, double agemax, double dateprev1, double dateprev2, int mobilav, double bage, double fage, int firstpass, int lastpass, double p[], int cptcoveff){
+ /* back1, year, month, day of starting backprojection
agemin, agemax range of age
dateprev1 dateprev2 range of dates during which prevalence is computed
anback2 year of end of backprojection (same day and month as back1).
@@ -8442,7 +8589,7 @@ set ter svg size 640, 480\nunset log y\n
*/
int yearp, stepsize, hstepm, nhstepm, j, k, cptcod, i, h, i1, k4, nres=0;
double agec; /* generic age */
- double agelim, ppij, ppi, yp,yp1,yp2,jprojmean,mprojmean,anprojmean;
+ double agelim, ppij, ppi, yp,yp1,yp2; /* ,jintmean,mintmean,aintmean;*/
double *popeffectif,*popcount;
double ***p3mat;
/* double ***mobaverage; */
@@ -8485,21 +8632,21 @@ set ter svg size 640, 480\nunset log y\n
}
hstepm=hstepm/stepm;
- yp1=modf(dateintmean,&yp);/* extracts integral of datemean in yp and
- fractional in yp1 */
- anprojmean=yp;
- yp2=modf((yp1*12),&yp);
- mprojmean=yp;
- yp1=modf((yp2*30.5),&yp);
- jprojmean=yp;
- if(jprojmean==0) jprojmean=1;
- if(mprojmean==0) jprojmean=1;
+ /* yp1=modf(dateintmean,&yp);/\* extracts integral of datemean in yp and */
+ /* fractional in yp1 *\/ */
+ /* aintmean=yp; */
+ /* yp2=modf((yp1*12),&yp); */
+ /* mintmean=yp; */
+ /* yp1=modf((yp2*30.5),&yp); */
+ /* jintmean=yp; */
+ /* if(jintmean==0) jintmean=1; */
+ /* if(mintmean==0) jintmean=1; */
i1=pow(2,cptcoveff);
if (cptcovn < 1){i1=1;}
- fprintf(ficresfb,"# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jprojmean,mprojmean,anprojmean,dateintmean,dateprev1,dateprev2);
- printf("# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jprojmean,mprojmean,anprojmean,dateintmean,dateprev1,dateprev2);
+ fprintf(ficresfb,"# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jintmean,mintmean,aintmean,dateintmean,dateprev1,dateprev2);
+ printf("# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jintmean,mintmean,aintmean,dateintmean,dateprev1,dateprev2);
fprintf(ficresfb,"#****** Routine prevbackforecast **\n");
@@ -8524,10 +8671,10 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficresfb," b%d%d",i,j);
fprintf(ficresfb," b.%d",j);
}
- for (yearp=0; yearp>=(anback2-anback1);yearp -=stepsize) {
+ for (yearp=0; yearp>=(anbackf-anbackd);yearp -=stepsize) {
/* for (yearp=0; yearp<=(anproj2-anproj1);yearp +=stepsize) { */
fprintf(ficresfb,"\n");
- fprintf(ficresfb,"\n# Back Forecasting at date %.lf/%.lf/%.lf ",jback1,mback1,anback1+yearp);
+ fprintf(ficresfb,"\n# Back Forecasting at date %.lf/%.lf/%.lf ",jbackd,mbackd,anbackd+yearp);
/* printf("\n# Back Forecasting at date %.lf/%.lf/%.lf ",jback1,mback1,anback1+yearp); */
/* for (agec=bage; agec<=agemax-1; agec++){ /\* testing *\/ */
for (agec=bage; agec<=fage; agec++){ /* testing */
@@ -8550,7 +8697,7 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficresfb,"\n");
for(j=1;j<=cptcoveff;j++)
fprintf(ficresfb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
- fprintf(ficresfb,"%.f %.f ",anback1+yearp,agec-h*hstepm/YEARM*stepm);
+ fprintf(ficresfb,"%.f %.f ",anbackd+yearp,agec-h*hstepm/YEARM*stepm);
for(i=1; i<=nlstate+ndeath;i++) {
ppij=0.;ppi=0.;
for(j=1; j<=nlstate;j++) {
@@ -8998,7 +9145,7 @@ void prwizard(int ncovmodel, int nlstate
/******************* Gompertz Likelihood ******************************/
double gompertz(double x[])
{
- double A,B,L=0.0,sump=0.,num=0.;
+ double A=0.0,B=0.,L=0.0,sump=0.,num=0.;
int i,n=0; /* n is the size of the sample */
for (i=1;i<=imx ; i++) {
@@ -9006,28 +9153,34 @@ double gompertz(double x[])
/* sump=sump+1;*/
num=num+1;
}
-
-
+ L=0.0;
+ /* agegomp=AGEGOMP; */
/* for (i=0; i<=imx; i++)
if (wav[i]>0) printf("i=%d ageex=%lf agecens=%lf agedc=%lf cens=%d %d\n" ,i,ageexmed[i],agecens[i],agedc[i],cens[i],wav[i]);*/
- for (i=1;i<=imx ; i++)
- {
- if (cens[i] == 1 && wav[i]>1)
- A=-x[1]/(x[2])*(exp(x[2]*(agecens[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp)));
-
- if (cens[i] == 0 && wav[i]>1)
+ for (i=1;i<=imx ; i++) {
+ /* mu(a)=mu(agecomp)*exp(teta*(age-agegomp))
+ mu(a)=x[1]*exp(x[2]*(age-agegomp)); x[1] and x[2] are per year.
+ * L= Product mu(agedeces)exp(-\int_ageexam^agedc mu(u) du ) for a death between agedc (in month)
+ * and agedc +1 month, cens[i]=0: log(x[1]/YEARM)
+ * +
+ * exp(-\int_ageexam^agecens mu(u) du ) when censored, cens[i]=1
+ */
+ if (wav[i] > 1 || agedc[i] < AGESUP) {
+ if (cens[i] == 1){
+ A=-x[1]/(x[2])*(exp(x[2]*(agecens[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp)));
+ } else if (cens[i] == 0){
A=-x[1]/(x[2])*(exp(x[2]*(agedc[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp)))
- +log(x[1]/YEARM)+x[2]*(agedc[i]-agegomp)+log(YEARM);
-
+ +log(x[1]/YEARM) +x[2]*(agedc[i]-agegomp)+log(YEARM);
+ } else
+ printf("Gompertz cens[%d] neither 1 nor 0\n",i);
/*if (wav[i] > 1 && agecens[i] > 15) {*/ /* ??? */
- if (wav[i] > 1 ) { /* ??? */
- L=L+A*weight[i];
+ L=L+A*weight[i];
/* printf("\ni=%d A=%f L=%lf x[1]=%lf x[2]=%lf ageex=%lf agecens=%lf cens=%d agedc=%lf weight=%lf\n",i,A,L,x[1],x[2],ageexmed[i]*12,agecens[i]*12,cens[i],agedc[i]*12,weight[i]);*/
- }
- }
+ }
+ }
- /*printf("x1=%2.9f x2=%2.9f x3=%2.9f L=%f\n",x[1],x[2],x[3],L);*/
+ /*printf("x1=%2.9f x2=%2.9f x3=%2.9f L=%f\n",x[1],x[2],x[3],L);*/
return -2*L*num/sump;
}
@@ -9036,7 +9189,7 @@ double gompertz(double x[])
/******************* Gompertz_f Likelihood ******************************/
double gompertz_f(const gsl_vector *v, void *params)
{
- double A,B,LL=0.0,sump=0.,num=0.;
+ double A=0.,B=0.,LL=0.0,sump=0.,num=0.;
double *x= (double *) v->data;
int i,n=0; /* n is the size of the sample */
@@ -9129,6 +9282,7 @@ int readdata(char datafile[], int firsto
int i=0, j=0, n=0, iv=0, v;
int lstra;
int linei, month, year,iout;
+ int noffset=0; /* This is the offset if BOM data file */
char line[MAXLINE], linetmp[MAXLINE];
char stra[MAXLINE], strb[MAXLINE];
char *stratrunc;
@@ -9162,8 +9316,53 @@ int readdata(char datafile[], int firsto
fprintf(ficlog,"Problem while opening datafile: %s with errno='%s'\n", datafile,strerror(errno));fflush(ficlog);return 1;
}
- i=1;
+ /* Is it a BOM UTF-8 Windows file? */
+ /* First data line */
linei=0;
+ while(fgets(line, MAXLINE, fic)) {
+ noffset=0;
+ if( line[0] == (char)0xEF && line[1] == (char)0xBB) /* EF BB BF */
+ {
+ noffset=noffset+3;
+ printf("# Data file '%s' is an UTF8 BOM file, please convert to UTF8 or ascii file and rerun.\n",datafile);fflush(stdout);
+ fprintf(ficlog,"# Data file '%s' is an UTF8 BOM file, please convert to UTF8 or ascii file and rerun.\n",datafile);
+ fflush(ficlog); return 1;
+ }
+ /* else if( line[0] == (char)0xFE && line[1] == (char)0xFF)*/
+ else if( line[0] == (char)0xFF && line[1] == (char)0xFE)
+ {
+ noffset=noffset+2;
+ printf("# Error Data file '%s' is a huge UTF16BE BOM file, please convert to UTF8 or ascii file (for example with dos2unix) and rerun.\n",datafile);fflush(stdout);
+ fprintf(ficlog,"# Error Data file '%s' is a huge UTF16BE BOM file, please convert to UTF8 or ascii file (for example with dos2unix) and rerun.\n",datafile);
+ fflush(ficlog); return 1;
+ }
+ else if( line[0] == 0 && line[1] == 0)
+ {
+ if( line[2] == (char)0xFE && line[3] == (char)0xFF){
+ noffset=noffset+4;
+ printf("# Error Data file '%s' is a huge UTF16BE BOM file, please convert to UTF8 or ascii file (for example with dos2unix) and rerun.\n",datafile);fflush(stdout);
+ fprintf(ficlog,"# Error Data file '%s' is a huge UTF16BE BOM file, please convert to UTF8 or ascii file (for example with dos2unix) and rerun.\n",datafile);
+ fflush(ficlog); return 1;
+ }
+ } else{
+ ;/*printf(" Not a BOM file\n");*/
+ }
+ /* If line starts with a # it is a comment */
+ if (line[noffset] == '#') {
+ linei=linei+1;
+ break;
+ }else{
+ break;
+ }
+ }
+ fclose(fic);
+ if((fic=fopen(datafile,"r"))==NULL) {
+ printf("Problem while opening datafile: %s with errno='%s'\n", datafile,strerror(errno));fflush(stdout);
+ fprintf(ficlog,"Problem while opening datafile: %s with errno='%s'\n", datafile,strerror(errno));fflush(ficlog);return 1;
+ }
+ /* Not a Bom file */
+
+ i=1;
while ((fgets(line, MAXLINE, fic) != NULL) &&((i >= firstobs) && (i <=lastobs))) {
linei=linei+1;
for(j=strlen(line); j>=0;j--){ /* Untabifies line */
@@ -9284,7 +9483,11 @@ int readdata(char datafile[], int firsto
return 1;
}
anint[j][i]= (double) year;
- mint[j][i]= (double)month;
+ mint[j][i]= (double)month;
+ /* if( (int)anint[j][i]+ (int)(mint[j][i])/12. < (int) (moisnais[i]/12.+annais[i])){ */
+ /* printf("Warning reading data around '%s' at line number %d for individual %d, '%s'\nThe date of interview (%2d/%4d) at wave %d occurred before the date of birth (%2d/%4d).\n",strb, linei,i, line, mint[j][i],anint[j][i], moisnais[i],annais[i]); */
+ /* fprintf(ficlog,"Warning reading data around '%s' at line number %d for individual %d, '%s'\nThe date of interview (%2d/%4d) at wave %d occurred before the date of birth (%2d/%4d).\n",strb, linei,i, line, mint[j][i],anint[j][i], moisnais[i],annais[i]); */
+ /* } */
strcpy(line,stra);
} /* End loop on waves */
@@ -9323,7 +9526,14 @@ int readdata(char datafile[], int firsto
}
annais[i]=(double)(year);
- moisnais[i]=(double)(month);
+ moisnais[i]=(double)(month);
+ for (j=1;j<=maxwav;j++){
+ if( (int)anint[j][i]+ (int)(mint[j][i])/12. < (int) (moisnais[i]/12.+annais[i])){
+ printf("Warning reading data around '%s' at line number %d for individual %d, '%s'\nThe date of interview (%2d/%4d) at wave %d occurred before the date of birth (%2d/%4d).\n",strb, linei,i, line, (int)mint[j][i],(int)anint[j][i], j,(int)moisnais[i],(int)annais[i]);
+ fprintf(ficlog,"Warning reading data around '%s' at line number %d for individual %d, '%s'\nThe date of interview (%2d/%4d) at wave %d occurred before the date of birth (%2d/%4d).\n",strb, linei,i, line, (int)mint[j][i],(int)anint[j][i], j, (int)moisnais[i],(int)annais[i]);
+ }
+ }
+
strcpy(line,stra);
/* Sample weight */
@@ -9343,6 +9553,8 @@ int readdata(char datafile[], int firsto
cutv(stra, strb, line, ' ');
if(strb[0]=='.') { /* Missing value */
lval=-1;
+ coqvar[iv][i]=NAN;
+ covar[ncovcol+iv][i]=NAN; /* including qvar in standard covar for performance reasons */
}else{
errno=0;
/* what_kind_of_number(strb); */
@@ -9446,7 +9658,6 @@ int decoderesult ( char resultline[], in
char stra[80], strb[80], strc[80], strd[80],stre[80];
removefirstspace(&resultline);
- printf("decoderesult:%s\n",resultline);
if (strstr(resultline,"v") !=0){
printf("Error. 'v' must be in upper case 'V' result: %s ",resultline);
@@ -9461,15 +9672,14 @@ int decoderesult ( char resultline[], in
TKresult[nres]=0; /* Combination for the nresult and the model */
return (0);
}
-
if( j != cptcovs ){ /* Be careful if a variable is in a product but not single */
- printf("ERROR: the number of variable in the resultline, %d, differs from the number of variable used in the model line, %d.\n",j, cptcovs);
- fprintf(ficlog,"ERROR: the number of variable in the resultline, %d, differs from the number of variable used in the model line, %d.\n",j, cptcovs);
+ printf("ERROR: the number of variables in the resultline, %d, differs from the number of variables used in the model line, %d.\n",j, cptcovs);
+ fprintf(ficlog,"ERROR: the number of variables in the resultline, %d, differs from the number of variables used in the model line, %d.\n",j, cptcovs);
}
for(k=1; k<=j;k++){ /* Loop on any covariate of the result line */
if(nbocc(resultsav,'=') >1){
cutl(stra,strb,resultsav,' '); /* keeps in strb after the first ' '
- resultsav= V4=1 V5=25.1 V3=0 strb=V3=0 stra= V4=1 V5=25.1 */
+ resultsav= V4=1 V5=25.1 V3=0 stra= V5=25.1 V3=0 strb= V4=1 */
cutl(strc,strd,strb,'='); /* strb:V4=1 strc=1 strd=V4 */
}else
cutl(strc,strd,resultsav,'=');
@@ -9494,7 +9704,9 @@ int decoderesult ( char resultline[], in
}
}
if(match == 0){
- printf("Error in result line: %d value missing; result: %s, model=%s\n",k1, resultline, model);
+ printf("Error in result line: V%d is missing in result: %s according to model=%s\n",k1, resultline, model);
+ fprintf(ficlog,"Error in result line: V%d is missing in result: %s according to model=%s\n",k1, resultline, model);
+ return 1;
}
}
}
@@ -9511,8 +9723,12 @@ int decoderesult ( char resultline[], in
}
if(match == 0){
printf("Error in result line: %d value missing; result: %s, model=%s\n",k1, resultline, model);
+ fprintf(ficlog,"Error in result line: %d value missing; result: %s, model=%s\n",k1, resultline, model);
+ return 1;
}else if(match > 1){
printf("Error in result line: %d doubled; result: %s, model=%s\n",k2, resultline, model);
+ fprintf(ficlog,"Error in result line: %d doubled; result: %s, model=%s\n",k2, resultline, model);
+ return 1;
}
}
@@ -10700,8 +10916,8 @@ int hPijx(double *p, int bage, int fage)
/* for (agedeb=fage; agedeb>=bage; agedeb--){ /\* If stepm=6 months *\/ */
for (agedeb=bage; agedeb<=fage; agedeb++){ /* If stepm=6 months and estepm=24 (2 years) */
/* nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm); /\* Typically 20 years = 20*12/6=40 *\/ */
- nhstepm=(int) rint((agedeb-ageminl)*YEARM/stepm); /* Typically 20 years = 20*12/6=40 */
- nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 */
+ nhstepm=(int) rint((agedeb-ageminl)*YEARM/stepm+0.1)-1; /* Typically 20 years = 20*12/6=40 or 55*12/24=27.5-1.1=>27 */
+ nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 or 28*/
/* nhstepm=nhstepm*YEARM; aff par mois*/
@@ -10793,7 +11009,14 @@ int main(int argc, char *argv[])
int *tab;
int mobilavproj=0 , prevfcast=0 ; /* moving average of prev, If prevfcast=1 prevalence projection */
- int backcast=0; /* defined as global for mlikeli and mle*/
+ /* double anprojd, mprojd, jprojd; /\* For eventual projections *\/ */
+ /* double anprojf, mprojf, jprojf; */
+ /* double jintmean,mintmean,aintmean; */
+ int prvforecast = 0; /* Might be 1 (date of beginning of projection is a choice or 2 is the dateintmean */
+ int prvbackcast = 0; /* Might be 1 (date of beginning of projection is a choice or 2 is the dateintmean */
+ double yrfproj= 10.0; /* Number of years of forward projections */
+ double yrbproj= 10.0; /* Number of years of backward projections */
+ int prevbcast=0; /* defined as global for mlikeli and mle, replacing backcast */
int mobilav=0,popforecast=0;
int hstepm=0, nhstepm=0;
int agemortsup;
@@ -10818,8 +11041,9 @@ int main(int argc, char *argv[])
double *epj, vepp;
double dateprev1, dateprev2;
- double jproj1=1,mproj1=1,anproj1=2000,jproj2=1,mproj2=1,anproj2=2000, dateproj1=0, dateproj2=0;
- double jback1=1,mback1=1,anback1=2000,jback2=1,mback2=1,anback2=2000, dateback1=0, dateback2=0;
+ double jproj1=1,mproj1=1,anproj1=2000,jproj2=1,mproj2=1,anproj2=2000, dateproj1=0, dateproj2=0, dateprojd=0, dateprojf=0;
+ double jback1=1,mback1=1,anback1=2000,jback2=1,mback2=1,anback2=2000, dateback1=0, dateback2=0, datebackd=0, datebackf=0;
+
double **ximort;
char *alph[]={"a","a","b","c","d","e"}, str[4]="1234";
@@ -11014,7 +11238,8 @@ int main(int argc, char *argv[])
noffset=noffset+3;
printf("# File is an UTF8 Bom.\n"); // 0xBF
}
- else if( line[0] == (char)0xFE && line[1] == (char)0xFF)
+/* else if( line[0] == (char)0xFE && line[1] == (char)0xFF)*/
+ else if( line[0] == (char)0xFF && line[1] == (char)0xFE)
{
noffset=noffset+2;
printf("# File is an UTF16BE BOM file\n");
@@ -11102,8 +11327,8 @@ int main(int argc, char *argv[])
}
if((num_filled=sscanf(line,"model=1+age%[^.\n]", model)) !=EOF){
if (num_filled != 1){
- printf("ERROR %d: Model should be at minimum 'model=1+age' %s\n",num_filled, line);
- fprintf(ficlog,"ERROR %d: Model should be at minimum 'model=1+age' %s\n",num_filled, line);
+ printf("ERROR %d: Model should be at minimum 'model=1+age+' instead of '%s'\n",num_filled, line);
+ fprintf(ficlog,"ERROR %d: Model should be at minimum 'model=1+age+' instead of '%s'\n",num_filled, line);
model[0]='\0';
goto end;
}
@@ -11789,8 +12014,8 @@ Interval (in months) between two waves:
ximort[i][j]=(i == j ? 1.0 : 0.0);
}
- /*p[1]=0.0268; p[NDIM]=0.083;*/
- /*printf("%lf %lf", p[1], p[2]);*/
+ p[1]=0.0268; p[NDIM]=0.083;
+ /* printf("%lf %lf", p[1], p[2]); */
#ifdef GSL
@@ -11916,9 +12141,9 @@ Interval (in months) between two waves:
printf("%f [%f ; %f]\n",p[i],p[i]-2*sqrt(matcov[i][i]),p[i]+2*sqrt(matcov[i][i]));
fprintf(ficlog,"%f [%f ; %f]\n",p[i],p[i]-2*sqrt(matcov[i][i]),p[i]+2*sqrt(matcov[i][i]));
}
- lsurv=vector(1,AGESUP);
- lpop=vector(1,AGESUP);
- tpop=vector(1,AGESUP);
+ lsurv=vector(agegomp,AGESUP);
+ lpop=vector(agegomp,AGESUP);
+ tpop=vector(agegomp,AGESUP);
lsurv[agegomp]=100000;
for (k=agegomp;k<=AGESUP;k++) {
@@ -11965,9 +12190,9 @@ Please run with mle=-1 to get a correct
stepm, weightopt,\
model,imx,p,matcov,agemortsup);
- free_vector(lsurv,1,AGESUP);
- free_vector(lpop,1,AGESUP);
- free_vector(tpop,1,AGESUP);
+ free_vector(lsurv,agegomp,AGESUP);
+ free_vector(lpop,agegomp,AGESUP);
+ free_vector(tpop,agegomp,AGESUP);
free_matrix(ximort,1,NDIM,1,NDIM);
free_ivector(dcwave,firstobs,lastobs);
free_vector(agecens,firstobs,lastobs);
@@ -12168,6 +12393,7 @@ Please run with mle=-1 to get a correct
fputs(line,stdout);
fputs(line,ficparo);
fputs(line,ficlog);
+ fputs(line,ficres);
continue;
}else
break;
@@ -12213,6 +12439,7 @@ Please run with mle=-1 to get a correct
fputs(line,stdout);
fputs(line,ficparo);
fputs(line,ficlog);
+ fputs(line,ficres);
continue;
}else
break;
@@ -12238,6 +12465,7 @@ Please run with mle=-1 to get a correct
fputs(line,stdout);
fputs(line,ficparo);
fputs(line,ficlog);
+ fputs(line,ficres);
continue;
}else
break;
@@ -12260,95 +12488,114 @@ Please run with mle=-1 to get a correct
}
/* Results */
+ endishere=0;
nresult=0;
+ parameterline=0;
do{
if(!fgets(line, MAXLINE, ficpar)){
endishere=1;
- parameterline=14;
+ parameterline=15;
}else if (line[0] == '#') {
/* If line starts with a # it is a comment */
numlinepar++;
fputs(line,stdout);
fputs(line,ficparo);
fputs(line,ficlog);
+ fputs(line,ficres);
continue;
}else if(sscanf(line,"prevforecast=%[^\n]\n",modeltemp))
parameterline=11;
- else if(sscanf(line,"backcast=%[^\n]\n",modeltemp))
+ else if(sscanf(line,"prevbackcast=%[^\n]\n",modeltemp))
parameterline=12;
- else if(sscanf(line,"result:%[^\n]\n",modeltemp))
+ else if(sscanf(line,"result:%[^\n]\n",modeltemp)){
parameterline=13;
+ }
else{
parameterline=14;
}
- switch (parameterline){
+ switch (parameterline){ /* =0 only if only comments */
case 11:
- if((num_filled=sscanf(line,"prevforecast=%d starting-proj-date=%lf/%lf/%lf final-proj-date=%lf/%lf/%lf mobil_average=%d\n",&prevfcast,&jproj1,&mproj1,&anproj1,&jproj2,&mproj2,&anproj2,&mobilavproj)) !=EOF){
- if (num_filled != 8) {
- printf("Error: Not 8 (data)parameters in line but %d, for example:prevforecast=1 starting-proj-date=1/1/1990 final-proj-date=1/1/2000 mobil_average=0\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
- fprintf(ficlog,"Error: Not 8 (data)parameters in line but %d, for example:prevforecast=1 starting-proj-date=1/1/1990 final-proj-date=1/1/2000 mov_average=0\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
- goto end;
- }
- fprintf(ficparo,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
+ if((num_filled=sscanf(line,"prevforecast=%d starting-proj-date=%lf/%lf/%lf final-proj-date=%lf/%lf/%lf mobil_average=%d\n",&prevfcast,&jproj1,&mproj1,&anproj1,&jproj2,&mproj2,&anproj2,&mobilavproj)) !=EOF && (num_filled == 8)){
+ fprintf(ficparo,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
printf("prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
fprintf(ficlog,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
fprintf(ficres,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
/* day and month of proj2 are not used but only year anproj2.*/
dateproj1=anproj1+(mproj1-1)/12.+(jproj1-1)/365.;
dateproj2=anproj2+(mproj2-1)/12.+(jproj2-1)/365.;
-
+ prvforecast = 1;
+ }
+ else if((num_filled=sscanf(line,"prevforecast=%d yearsfproj=%lf mobil_average=%d\n",&prevfcast,&yrfproj,&mobilavproj)) !=EOF){/* && (num_filled == 3))*/
+ printf("prevforecast=%d yearsfproj=%.2lf mobil_average=%d\n",prevfcast,yrfproj,mobilavproj);
+ fprintf(ficlog,"prevforecast=%d yearsfproj=%.2lf mobil_average=%d\n",prevfcast,yrfproj,mobilavproj);
+ fprintf(ficres,"prevforecast=%d yearsfproj=%.2lf mobil_average=%d\n",prevfcast,yrfproj,mobilavproj);
+ prvforecast = 2;
+ }
+ else {
+ printf("Error: Not 8 (data)parameters in line but %d, for example:prevforecast=1 starting-proj-date=1/1/1990 final-proj-date=1/1/2000 mobil_average=0\nnor 3 (data)parameters, for example:prevforecast=1 yearsfproj=10 mobil_average=0. Your line=%s . You are running probably an older format.\n, ",num_filled,line);
+ fprintf(ficlog,"Error: Not 8 (data)parameters in line but %d, for example:prevforecast=1 starting-proj-date=1/1/1990 final-proj-date=1/1/2000 mobil_average=0\nnor 3 (data)parameters, for example:prevforecast=1 yearproj=10 mobil_average=0. Your line=%s . You are running probably an older format.\n, ",num_filled,line);
+ goto end;
}
break;
case 12:
- /*fscanf(ficpar,"backcast=%d starting-back-date=%lf/%lf/%lf final-back-date=%lf/%lf/%lf mobil_average=%d\n",&backcast,&jback1,&mback1,&anback1,&jback2,&mback2,&anback2,&mobilavproj);*/
- if((num_filled=sscanf(line,"backcast=%d starting-back-date=%lf/%lf/%lf final-back-date=%lf/%lf/%lf mobil_average=%d\n",&backcast,&jback1,&mback1,&anback1,&jback2,&mback2,&anback2,&mobilavproj)) !=EOF){
- if (num_filled != 8) {
- printf("Error: Not 8 (data)parameters in line but %d, for example:backcast=1 starting-back-date=1/1/1990 final-back-date=1/1/1970 mobil_average=1\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
- fprintf(ficlog,"Error: Not 8 (data)parameters in line but %d, for example:backcast=1 starting-back-date=1/1/1990 final-back-date=1/1/1970 mobil_average=1\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
- goto end;
- }
- printf("backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
- fprintf(ficparo,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
- fprintf(ficlog,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
- fprintf(ficres,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
- /* day and month of proj2 are not used but only year anproj2.*/
+ if((num_filled=sscanf(line,"prevbackcast=%d starting-back-date=%lf/%lf/%lf final-back-date=%lf/%lf/%lf mobil_average=%d\n",&prevbcast,&jback1,&mback1,&anback1,&jback2,&mback2,&anback2,&mobilavproj)) !=EOF && (num_filled == 8)){
+ fprintf(ficparo,"prevbackcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevbcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ printf("prevbackcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevbcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ fprintf(ficlog,"prevbackcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevbcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ fprintf(ficres,"prevbackcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevbcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ /* day and month of back2 are not used but only year anback2.*/
dateback1=anback1+(mback1-1)/12.+(jback1-1)/365.;
dateback2=anback2+(mback2-1)/12.+(jback2-1)/365.;
+ prvbackcast = 1;
+ }
+ else if((num_filled=sscanf(line,"prevbackcast=%d yearsbproj=%lf mobil_average=%d\n",&prevbcast,&yrbproj,&mobilavproj)) ==3){/* && (num_filled == 3))*/
+ printf("prevbackcast=%d yearsbproj=%.2lf mobil_average=%d\n",prevbcast,yrbproj,mobilavproj);
+ fprintf(ficlog,"prevbackcast=%d yearsbproj=%.2lf mobil_average=%d\n",prevbcast,yrbproj,mobilavproj);
+ fprintf(ficres,"prevbackcast=%d yearsbproj=%.2lf mobil_average=%d\n",prevbcast,yrbproj,mobilavproj);
+ prvbackcast = 2;
+ }
+ else {
+ printf("Error: Not 8 (data)parameters in line but %d, for example:prevbackcast=1 starting-back-date=1/1/1990 final-back-date=1/1/2000 mobil_average=0\nnor 3 (data)parameters, for example:prevbackcast=1 yearsbproj=10 mobil_average=0. Your line=%s . You are running probably an older format.\n, ",num_filled,line);
+ fprintf(ficlog,"Error: Not 8 (data)parameters in line but %d, for example:prevbackcast=1 starting-back-date=1/1/1990 final-back-date=1/1/2000 mobil_average=0\nnor 3 (data)parameters, for example:prevbackcast=1 yearbproj=10 mobil_average=0. Your line=%s . You are running probably an older format.\n, ",num_filled,line);
+ goto end;
}
break;
case 13:
- if((num_filled=sscanf(line,"result:%[^\n]\n",resultline)) !=EOF){
- if (num_filled == 0){
- resultline[0]='\0';
- printf("Warning %d: no result line! It should be at minimum 'result: V2=0 V1=1 or result:.\n%s\n", num_filled, line);
- fprintf(ficlog,"Warning %d: no result line! It should be at minimum 'result: V2=0 V1=1 or result:.\n%s\n", num_filled, line);
- break;
- } else if (num_filled != 1){
- printf("ERROR %d: result line! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",num_filled, line);
- fprintf(ficlog,"ERROR %d: result line! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",num_filled, line);
- }
- nresult++; /* Sum of resultlines */
- printf("Result %d: result=%s\n",nresult, resultline);
- if(nresult > MAXRESULTLINES){
- printf("ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
- fprintf(ficlog,"ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
- goto end;
- }
- decoderesult(resultline, nresult); /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */
+ num_filled=sscanf(line,"result:%[^\n]\n",resultline);
+ nresult++; /* Sum of resultlines */
+ printf("Result %d: result:%s\n",nresult, resultline);
+ if(nresult > MAXRESULTLINES){
+ printf("ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\nYou can use the 'r' parameter file '%s' which uses option mle=0 to get other results. ",MAXRESULTLINES,nresult,rfileres);
+ fprintf(ficlog,"ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\nYou can use the 'r' parameter file '%s' which uses option mle=0 to get other results. ",MAXRESULTLINES,nresult,rfileres);
+ goto end;
+ }
+ if(!decoderesult(resultline, nresult)){ /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */
fprintf(ficparo,"result: %s\n",resultline);
fprintf(ficres,"result: %s\n",resultline);
fprintf(ficlog,"result: %s\n",resultline);
- break;
- case 14:
- if(ncovmodel >2 && nresult==0 ){
- printf("ERROR: no result lines! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",line);
- goto end;
- }
- break;
- default:
- nresult=1;
- decoderesult(".",nresult ); /* No covariate */
+ } else
+ goto end;
+ break;
+ case 14:
+ printf("Error: Unknown command '%s'\n",line);
+ fprintf(ficlog,"Error: Unknown command '%s'\n",line);
+ if(line[0] == ' ' || line[0] == '\n'){
+ printf("It should not be an empty line '%s'\n",line);
+ fprintf(ficlog,"It should not be an empty line '%s'\n",line);
+ }
+ if(ncovmodel >=2 && nresult==0 ){
+ printf("ERROR: no result lines! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",line);
+ fprintf(ficlog,"ERROR: no result lines! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",line);
}
+ /* goto end; */
+ break;
+ case 15:
+ printf("End of resultlines.\n");
+ fprintf(ficlog,"End of resultlines.\n");
+ break;
+ default: /* parameterline =0 */
+ nresult=1;
+ decoderesult(".",nresult ); /* No covariate */
} /* End switch parameterline */
}while(endishere==0); /* End do */
@@ -12365,11 +12612,44 @@ This is probably because your parameter
Please run with mle=-1 to get a correct covariance matrix.\n",ageminpar,agemaxpar);
}else{
/* printinggnuplot(fileresu, optionfilefiname,ageminpar,agemaxpar,fage, prevfcast, backcast, pathc,p, (int)anproj1-(int)agemin, (int)anback1-(int)agemax+1); */
- printinggnuplot(fileresu, optionfilefiname,ageminpar,agemaxpar,bage, fage, prevfcast, backcast, pathc,p, (int)anproj1-bage, (int)anback1-fage);
+ /* It seems that anprojd which is computed from the mean year at interview which is known yet because of freqsummary */
+ /* date2dmy(dateintmean,&jintmean,&mintmean,&aintmean); */ /* Done in freqsummary */
+ if(prvforecast==1){
+ dateprojd=(jproj1+12*mproj1+365*anproj1)/365;
+ jprojd=jproj1;
+ mprojd=mproj1;
+ anprojd=anproj1;
+ dateprojf=(jproj2+12*mproj2+365*anproj2)/365;
+ jprojf=jproj2;
+ mprojf=mproj2;
+ anprojf=anproj2;
+ } else if(prvforecast == 2){
+ dateprojd=dateintmean;
+ date2dmy(dateprojd,&jprojd, &mprojd, &anprojd);
+ dateprojf=dateintmean+yrfproj;
+ date2dmy(dateprojf,&jprojf, &mprojf, &anprojf);
+ }
+ if(prvbackcast==1){
+ datebackd=(jback1+12*mback1+365*anback1)/365;
+ jbackd=jback1;
+ mbackd=mback1;
+ anbackd=anback1;
+ datebackf=(jback2+12*mback2+365*anback2)/365;
+ jbackf=jback2;
+ mbackf=mback2;
+ anbackf=anback2;
+ } else if(prvbackcast == 2){
+ datebackd=dateintmean;
+ date2dmy(datebackd,&jbackd, &mbackd, &anbackd);
+ datebackf=dateintmean-yrbproj;
+ date2dmy(datebackf,&jbackf, &mbackf, &anbackf);
+ }
+
+ printinggnuplot(fileresu, optionfilefiname,ageminpar,agemaxpar,bage, fage, prevfcast, prevbcast, pathc,p, (int)anprojd-bage, (int)anbackd-fage);
}
printinghtml(fileresu,title,datafile, firstpass, lastpass, stepm, weightopt, \
- model,imx,jmin,jmax,jmean,rfileres,popforecast,mobilav,prevfcast,mobilavproj,backcast, estepm, \
- jprev1,mprev1,anprev1,dateprev1, dateproj1, dateback1,jprev2,mprev2,anprev2,dateprev2,dateproj2, dateback2);
+ model,imx,jmin,jmax,jmean,rfileres,popforecast,mobilav,prevfcast,mobilavproj,prevbcast, estepm, \
+ jprev1,mprev1,anprev1,dateprev1, dateprojd, datebackd,jprev2,mprev2,anprev2,dateprev2,dateprojf, datebackf);
/*------------ free_vector -------------*/
/* chdir(path); */
@@ -12442,13 +12722,23 @@ Please run with mle=-1 to get a correct
}/* end if moving average */
/*---------- Forecasting ------------------*/
- if(prevfcast==1){
- /* if(stepm ==1){*/
- prevforecast(fileresu, anproj1, mproj1, jproj1, agemin, agemax, dateprev1, dateprev2, mobilavproj, mobaverage, bage, fage, firstpass, lastpass, anproj2, p, cptcoveff);
+ if(prevfcast==1){
+ /* /\* if(stepm ==1){*\/ */
+ /* /\* anproj1, mproj1, jproj1 either read explicitly or yrfproj *\/ */
+ /*This done previously after freqsummary.*/
+ /* dateprojd=(jproj1+12*mproj1+365*anproj1)/365; */
+ /* dateprojf=(jproj2+12*mproj2+365*anproj2)/365; */
+
+ /* } else if (prvforecast==2){ */
+ /* /\* if(stepm ==1){*\/ */
+ /* /\* anproj1, mproj1, jproj1 either read explicitly or yrfproj *\/ */
+ /* } */
+ /*prevforecast(fileresu, dateintmean, anproj1, mproj1, jproj1, agemin, agemax, dateprev1, dateprev2, mobilavproj, mobaverage, bage, fage, firstpass, lastpass, anproj2, p, cptcoveff);*/
+ prevforecast(fileresu,dateintmean, dateprojd, dateprojf, agemin, agemax, dateprev1, dateprev2, mobilavproj, mobaverage, bage, fage, firstpass, lastpass, p, cptcoveff);
}
- /* Backcasting */
- if(backcast==1){
+ /* Prevbcasting */
+ if(prevbcast==1){
ddnewms=matrix(1,nlstate+ndeath,1,nlstate+ndeath);
ddoldms=matrix(1,nlstate+ndeath,1,nlstate+ndeath);
ddsavms=matrix(1,nlstate+ndeath,1,nlstate+ndeath);
@@ -12463,8 +12753,14 @@ Please run with mle=-1 to get a correct
hBijx(p, bage, fage, mobaverage);
fclose(ficrespijb);
- prevbackforecast(fileresu, mobaverage, anback1, mback1, jback1, agemin, agemax, dateprev1, dateprev2,
- mobilavproj, bage, fage, firstpass, lastpass, anback2, p, cptcoveff);
+ /* /\* prevbackforecast(fileresu, mobaverage, anback1, mback1, jback1, agemin, agemax, dateprev1, dateprev2, *\/ */
+ /* /\* mobilavproj, bage, fage, firstpass, lastpass, anback2, p, cptcoveff); *\/ */
+ /* prevbackforecast(fileresu, mobaverage, anback1, mback1, jback1, agemin, agemax, dateprev1, dateprev2, */
+ /* mobilavproj, bage, fage, firstpass, lastpass, anback2, p, cptcoveff); */
+ prevbackforecast(fileresu, mobaverage, dateintmean, dateprojd, dateprojf, agemin, agemax, dateprev1, dateprev2,
+ mobilavproj, bage, fage, firstpass, lastpass, p, cptcoveff);
+
+
varbprlim(fileresu, nresult, mobaverage, mobilavproj, bage, fage, bprlim, &ncvyear, ftolpl, p, matcov, delti, stepm, cptcoveff);
@@ -12472,7 +12768,7 @@ Please run with mle=-1 to get a correct
free_matrix(ddnewms, 1, nlstate+ndeath, 1, nlstate+ndeath);
free_matrix(ddsavms, 1, nlstate+ndeath, 1, nlstate+ndeath);
free_matrix(ddoldms, 1, nlstate+ndeath, 1, nlstate+ndeath);
- } /* end Backcasting */
+ } /* end Prevbcasting */
/* ------ Other prevalence ratios------------ */